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OH AP 'HER I 


EilROLUCTION 

1.1 FOUR AREAS OE MATHEMATICAL BIOSCIENCES CONSIDERED IN 
THE PRESENT THESIS 

In the present thesis, we have made contribution to 
the following four areas of Mathematical Biosciences, namely, 

(i) optimal exploitation of renewable resources like 
fisheries 

(ii) optimal exploitation of exhaustible resources like mines 

(iii) linear birth-death-immigration-emigration processes 

(iv) spread of epidemics. 

In each of these areas, we have investigated the 
solution of a number of problems. Thus in optimal exploitation 
of fisheries, we consider 

(i) use of dynamic programming, calculus of variations and 
maximum principle to find the maximum value of the 
present value of total profits over an infinite span 
of time for a single species 

(ii) use of same techniques for the optimal exploitation of 
multispecies communities 

(iii) obtaining of limiting optimal solutions when overall 
optimal solutions are not obtainable. 

In optimal exploitation of mines, we consider two 


problems 
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(i) optimal exploitation when time span is fixed 

(ii) optimal exploitation when time span is variable and 
can be chosen so as to maximize the maximum present 
value of the profit for variation in length of time span. 

In linear birth-death- immigration-emigration processes, 

we consider 

(i) linear birth-death-emigration process (unsteady case) 

(ii) discrete birth-death-immigrat ion-emigration process . 

In epidemic theory, we consider the spread of epidemic 
for the case of multiple infections and study 

(i) deterministic model 

(ii) stochastic model 

Each of these areas has a large literature including 
dozens of standard level monographs and hundreds of research 
papers. In the next four sections we give a brief survey of 
that part of the literature which is related to our work in 
the present thesis in order to give a perspective for the work 
we have done. We also give a very brief summary of the results 
obtained in all areas. Chapter II - V will give detailed 
derivation of the results. 

Most of the work reported in the thesis has been 
published or accepted for publication or is ready for communica- 
tion in the form of the following papers s 

Optimal exploitation of controlled fishers, to appear 
in-J. Math. Physi. Sci. 
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Optimal exploitation of multispecies fisheries, published, 
in Some Aspects of Mechanics of Gontinua, part 1, 1977. 
Limiting optimal solutions for exploitation of fisheries, 
to appear in Proc. Nat. Acad. Sci. 

A Mathematical model for optimal exploitation of mines, 
to appear in J . Math. Physi. Sci. 

A class of birth— death-emigration processes, pub. in 
Nat, Acad. Sci. letters, Yol. 2, 1, 1979. 

Discrete birth-death-emigration-immigration processes, 
to be communicated. 

1.2. OPTIMAL EXPLOITATION OP RENEWABLE RESOURCES : PISHEEUES 
1.2.1 SURYEY OP LITERATURE 

Renewable resources are those capable of regeneration 
like fish, timber etc. Clark £ 1976]] states : " they grow by 
the gift of nature.” 

The most important feature of a renewable resource also 
called replenishable or reproducible resource is its law of 
biological growth. This characteristic can be studied in many 
ways e.g. as growth in the number of members vide Berverton and 
Hott 1957]] , or as growth in the total biomass of the species 
vide Lotka £ 1956]]. 

Extremely simple mathematical models for the exploita- 
tion of biological resources (renewable resources) are given by 
the. ordinary differential equation of the form £ Clark [1 1976]] ]] 



where 


x = x(t) denotes the size of the population or 
biomass at time t. 

f(x) = a given function defining the natural growth 
rate of p opulat ion * 

and h(t) = rate of removal or harvesting. 

The management of renewable resources had been based 
on the concept of maximum sustainable yield. This concept was 
based on the assumption that at any population level, less than 
a certain level, a surplus of production existed which could be 
harvested without disturbing the stock level. The models 
exhibiting this trend of growth were studied by G-ordom [[ 1954]], 
Shaefer £.1957]], Plourde [][l97o][] and Christy [[ 1973]]]]. But 
this concept had to face several objections raised on both 
biological and socio-economic grounds. Different forms of 
growth functions have been reported by May [[ 1973 ]] and 
Clark C 1976] . 

The subject of an economic theory, developed in the 
literature since Gordon [[ 1953]], has been the economic effects 
of uncontrolled exploitation of a common property resource. 
Gordon [[ 1954:]]] and Scott [1 1955]] have made contributions to 
the common property resources. Such resources where ownership 
is not perfect, are often over exploited since they allow the 
entry of a surplus of labour and capital into the industry. 



5 


In his articles Gordon Cl9£3,54 3 argued that this entry would 
continue untill all rent of fishery was dissipated wide 
Clark £1976]]. A theory exists in support of the argument 
that this type of exploitation would lead to the dissipation 
of profits and ultimately depletion of the resources. 

Thus for policy makers the problem has been to device 
the methods to correct this misallocat ion of resources. Apart 
from other ways, the greater control of exploitation has been 
done by putting restrictions on quotas and efforts. Remaining 
in the domain of quota or effort restrictions, fisheries 
exploitation can be managed to attain various objectives which 
may be economic, social or both. 

A number of contributions to the economic analysis of 
the exploitation of fisheries have been due to many including 
Smith [1 1968]], Scott £ 1955 ]] , Clark £ 1971, 72, 73 £ , Bedding ton, 
Watts and Wright £ 1975]], Waugh and Calvo £ 1974]] and 
Colbert and Swan £ 1975]]. Smith £ 1968]] had provided a 
unified general model which with simplifying assumptions could 
be applicable to a wide range of resource problems. This 
general model could not make out the detailed investigations 
of question of extinction as well as levels and stability of 
equilibriums. Particular questions concerning the behaviour 
of the population as it nears zero or critical levels of 
exploitation were considered in detail by Clark £ 1971, 72]], 
Clark £ 1973 J studied a mathematical model for the Antartic 
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■blue whale with constant price and a particular function for 
total costs. Colbert and Swan £ 19753 made good replacement 
of weak assumptions of Clark Q 1973 3 regarding price and cost 
function. 

Clark Q 19 71, 72 3 had taken constant costs as well as 
given prices which prevented him from assessing the effects 
of price elasticities and of different characteristics of 
production function. Beddington, Walts and W right Q 1975 3 
had studied the consequences of different forms of production 
functions. The object was to maximize the present value of 
total profits subject to a production function and two differen- 
tial equations defining the rate of change of capital stock and 
population. Waugh and Calvo £ 1974]] saw the fisheries as 
exhaustible resources with the same object. 

Most of the contributions to the economic theory of 
exploitation of fisheries has been through the maximization 
of some objective function which is economic in nature. A 
very commonly used objective function is the present value of 
total profits over time. This is maximized subject to a 
differential equation representing the rate of change of 
population along with a restriction on quota. 

Sancho and Mitchell £ 1975 3 used dynamic programming 
technique to obtain the optimal effort each year to optimize 
the present value of total profits over time, Sancho and 
Mitchell [1 1977^] applied the model to the Canadian ground fish 
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industry and using the data over ten-year period, obtained the 
optimal effort as well as optimal present value of total profits. 
Silvert and Smith Q 1977 ^ considered the maximization of 
discounted present value of net Revenue given by a function of 
vector population and vector harvest while studying the 
exploitation of multispecies community* 

1.2,2 DESCRIP II OR OR PROBIEMS MB RESUMS 

Sancho and Mitchell |~ 19 75 ^ used dynamic programming 
to obtain optimal effort each year to optimize the present 
value of total profits over time. Subject to a differential 
equation defining the rate of change of fish population and a 
constraint on effort in order to take into account the specified 
quota that may be alloted to the fishing party or the capacity 
of the market to buy fish. Sancho and Mitchell £ 1975]] found 
their results in terms of R, the initial p op ulat ion. Thus if 
the maximum value of total profits is needed each year, we have 
to make the estimate of R each year which can be done heuris- 
tically by taking the ratio of catch to the effort and comparing 
one with the next. 

In our first problem of optimal exploitation of 
controlled fisheries, firstly we study the model of Sancho and 
Mitchell £1975]] by using Maximum principle and Calculus of 
variations. Then we impose the constraint on catch per unit 
time. Finally, considering exploitation over finite horizon of 
time the following results are found out. 



We are able to obtain our results explicitly as 
functions of time and thus avoid the heuristic 
approach, which is necessary in Sancho and Mitchell 
L 1975 3 if the maximum value of total profits is 
needed each year. 

We can also find out the conditions under which the 
solutions are optimal. 

We can also find out the conditions which are to be 
satisfied if we impose the constraint on catch per 
unit t ime . 

Exploitation over finite horizon of time does not 
requires the above conditions. 

In the second problem, we use optimality principle 
of dynamic programming and the maximum principle to optimize 
the present value of total profits (P.V) over all time in a 
fishery exploitation when the fishery has a number of different 
but non-interacting n species. This is because the area of 
the pond or a river or a sea may contain a number of different 
types of fish commanding different prices and all of them may 
be caught. We get the following results. 

The maximum value of pre'sent value of profits is a 
quadratic function of R^,Rg, . . . ,R n , the initial fish 
populations. As an example, a particular case of two 
species is solved completely. 

The case when the catch of any species per unit time 
depends only on its own population is solved by using 
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maximum principle still getting present value as 
a quadratic function of R^,Hg, . . .,R n . Some more 
results for the limiting case are also pointed out. 

In the last problem, it is shown that limiting optimal 
solutions may he possible for an ecological system which may 
not, otherwise, be tractable. We consider a general system. 

As a special case the system considered by Clark L_ 1976H is 
solved completely for optimal equilibrium solutions. She case 
for interacting species is also considered. We observe 
When we are not able to find the optimal non- 
equilibrium singular controls, we find the equilibrium 
optimal solutions. If initial populations are less 
than the equilibrium values we don’t harvest till 
populations rise to equilibrium values. If either or 
both initial populations are more than equilibrium 
values, we do maximum harvesting possible required 
by the market till equilibrium values are reached. 
After equilibrium values are reached, we use the 
optimal effort (limiting value) thereafter and the 
populations remain stable. 

1.3 QPTB1AL EXPLOITATION OF EXHAUSTIBLE RESOURCES ; MINES 
1.3 .i SURVEY OF LITERATURE 

Exhaustible resources are those not capable of 


regeneration, such as petroleum, natural gas and minerals etc. 
They are fixed assets. As industry has developed, the need for 
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minerals lias increased and the search for more minerals has 
gone on. Today, the mineral resources of the earth are so 
commonly used that they affect every aspect of our lives. 

Minerals have some distinct characteristics which 
distinguish them from other natural resources. According to 
Erich Zimmerman, ” the most disconcerting feature of minerals 
is their exhaust ibility, " By and large, minerals are depleting 
assets and once mined out the deposits get depleted without any 
chance of replenishment. 

The mathematical treatment of the subject goes back 
to Hotelling [[3 1931 ^ wll ° emphasized the necessity of employing 
the methods of calculus of variations in the theory of exhaus- 
tible resources. The earlier works by Gray |~ 19143 and 
Hotelling Q 19313 concerning the M The theory of mine" have 
been extended by many including Smith |~1968^, Waugh and 
Oalvo [1 1974 3, Gordon £ 1967 3, Goldsmith £ 1974 3 and 
Glark j3 1976 3 • 

But in all these studies, while dealing with exhaus- 
tible resources, the effort that we apply while extracting the 
amount of mineral per unit time has not been taken into account. 
This aspect is studied in the third chapter by assuming different 
extraction functions and with some cost function as also 
functions of effort applied per unit time. 

1.3.2 DESCRIPTION OE PROBLEMS AND RESULTS 

In all the problems, our aim is to maximize the 
present value of total profits over time (P.Y). The main 
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difference from the renewable resources is that the mine is to 
be exhausted in a finite time- interval. We deal with two types 
of problems, in the first we maximize present value of total 
profits so that the mine is exhausted in a given time- interval 
and in the second we choose that t ime- interval for exploitation 
which provides us with the maximum of all the maximum P.Vs. 
'fating cost function Quadratic in effort QSancho and Mitchell 
|” 19 75 21 lithe procedure is illustrated numerically in all the 
problems which are given below along with results. 

firstly, a simple mathematical model is used in which 
the extraction function is linear in amount of material and the 
effort. Calculus of variations and Maximum principle are used 
to obtain the amount of the material to be taken out and the 
optimal effort for getting the maximum of the P.7, over the 
time-period in which the mine is exhausted. We observe the 
following 

for a fixed set of values of parameters the maximum 
profit increases with time and finally stabilizes 
showing that in very large time we do not get any 
additional profit. Graphs for amount of material and 
optimal effort against time are also plotted. 

Secondly, we consider the extraction function as 
Quadratic in amount of material and as linear in effort. Using 
calculus of variations, we deal w it h the same problem as above. 
Giving a numerical example, the following result is observed. 
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Slower exploitation of a mine over a long period 
is more profitable than quicker exploitation over a 
short time. The variations of amount of material, 
the optimal effort and the rate of change of the 
amount of material with time are shown in different 
graphs. 

lastly, we consider the problem in which the extraction 
function is some product of the amount of material in the mine 
and the effort. Using calculus of variations we are again 
faced with the similar problem as above. The result is as 
follows : 

Even in this case slow exploitation is preferred to 
quicker exploitation since the latter is more profi- 
table. Graphs show the variation of the amount of 
material and the effort with time, 

1.4 LINEAR BIRTH-DEATH- IMM IGRAIIOU— EMIGRATION PROCESSES 

1.4.1 SURVEY OP LITERATURE 

Early studies in the mathematical theory of population 
growth, epidemics, queues and several other situations were 
primarily 'concerned with the development of deterministic models 
which determine the number of individuals in the system at any 
time t. In this approach, on the basis of some assumptions 
concerning changes in numbers of individuals in the system in 
a small interval of time, an equation (differential or integral, 
etc.) for the number of individuals in the system, say x(t), is 
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formed. This equation is solved subject to some initial 
condition to find out x(t). By assuming x(t) as a real valued 
continuous function of time and the change in the small interval 
of time as deterministic, one obtains the result that number of 
individuals in the system does not change if initial conditions 
are not altered. Therefore, this approach does not take into 
consideration the large number of random or chance factors that 
can influence the system. 

In 1939, birth and death processes were introduced by 
W. feller. These processes have since been used to describe 
stochastic models for population growth, queues, epidemics 
and telephone calls etc. la this approach, a random variable 
X(t), integer (or real) valued, is used to describe the size 
of the system at time t and the change in the system in a 
small interval of time is expressed in terms of probabilities 
of certain elementary events in that interval, let P(n,t) 
denote the probability of there being n persons in the system 
at time t i.e. 

P(n,t) = PQx(t) = n^J , n = 0,1, 2,3,... 

In birth-and-death processes, a system of differential-difference 
equations in P(n,t) is formulated. Then it is solved with some 
initial conditions to find out the size of the system at any 
t ime t. 

When the population size is n, let the probability of 
a birth, a death, an immigration and an emigration in time 
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interval (t, t + At) be given respectively by A At + O(At), 

A t + 0( At ) , v n At + 0( At) and a n At + O(At). 

Depending upon the nature of the probability function, 
the following processes have been classified in the literature, 

( i) Poisson Process when A = A , u = 0, v = 0 and a =0 

n n n n 

(ii) Pure Birth Process «’ A^ = A q , n n = 0, v n = o and a n = 0 

(iii) Pure Death Process A = 0, u = y , v =0 and a =0 

30 30 30, 30 30 

(iv) Linear Birth-Death 

Process ” A n = nA , w n = ny , v n = 0 and a n = 0 

(v) Don-Homogenous 

Birth-Death Process " A =A (t),u =u (t),v ==0 and a =0 

O iLX OX XX XX XX 

(vi) Linear Birth-Death- 

Immigrat ion Process M A n =nA,u n =nv, v n = v and a n = 0 

(vii) Linear Birth-Death— 

Immigr at ion-Em igr at ion 

Process ** A =nA,y =ny , v = v and a =* a. 

ii IX JO 30 

'The complete solutions for the first six cases (when 
a n = 0) have been given, among others, by Bailey [1 1964 3? 
Bartlett [[ 1955 3 > Bharucha-Re id 3 I960 3 , Peller j3 1959 3 ? 
Iosifescu and Tautu j3 19733? Kendall £ 1948, 1952 3 and 
Srinivasan and Mehta [3 1976 3 . 

Recently Getz [[ 1975,19763 attempted to discuss the 
case (vii) when all the four processes are involved, but he was 
led to a deceptively simple treatment by omitting a term of 
probability of extinction. This error was pointed out by 
Kapur [[I979a3. Kapur [[ 1979b, 1979c3 also found out the 
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complete solution in terms of incomplete beta functions when 
three processes i.e. birth, death and emigration are involved. 
Kapur jjLS78a, 78b, ?8c, 78d,79d, 79£] and Kapur and Kapur jjL975]also 
obtained the complete solution for the steady case in terms 
of hypergeometr ic and generalized hypergeometr ic functions. 

1.4.2 DESCRIPTION OP PROBLEMS AND RESULTS 

Yfhile dealing with the birt3n.de at h-immigr at ion- 
emigration process Getz £ 1975, 1976^ found a partial differen- 
tial equation for generating function in which he omitted one 
term containing the probability of extinction. This omission 
makes all his calculations wrong. Knowing the probability of 
extinction is the foremost necessity if we want to know all 
about the system. Kapur £ 1979b, 1979c J calculated the 
probability of extinction of B-D-E process in terms of incomplete 
beta functions. We give here an alternative method for finding 
the probability of extinction. We put our problem as follows. 

The probability of extinction for a birth-death- 
emigration process for the case (a < u , emigration rate is less 
than the death rate) is obtained in terms of the solution of 
an Abel’s integral equation which is solved numerically. We 
observe 

As expected the agreement between our solution and 

that of Kapur [1 1979b, 19790^ is quite good. 

In another problem, we consider the discrete version 
of the continuous-time birth-death- immigration emigration process 
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and solve it numerically. The result is stated as follows 
We get more insight into the system and draw some 
general conclusions. Different probabilities have 
been plotted for n o = 3 and n Q = 6, where n Q is the 
initial number of persons in the system against time 
as well as against number of persons. Graphs for 
Mean and variance have also been plotted. 

1.5 SPREAD OR EPIDEMICS 

1.5.1 SURVEY OP LITERATURE 

Por the literature in the theory of epidemics, upto 
1957, an excellent monograph is by Bailey £ 1957]] . A 
comprehensive study of the mathematical contribution to the 
description of the spread of epidemics has been done by 
Dietz £19672] in the form of an excellent review. 

In recent years, much of the work in epidemics is 
related to the control problems studied by many including 
Bobisud £ 19772], Wickwire £ 19762] and Jaquette £l97o2] but 
we will not deal with these types of problems. 

Earlier, in the theory of epidemics, the disease was 
supposed to be spread by infected persons, called injectives, 
who were recognised by their having some symptoms. But Weiss 
£ 19652] introduced anew aspect in the spread of epidemics i.e. 
spread by carriers* A carrier is defined to be an individual 
who does not have overt disease symptom but who is however able 
to communicate the disease to others. Other contributions in 
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this field are due to Downton [_ 1968 J and Abakuks |_ 1974 ^ . 

Me will not also discuss this aspect . 

A widely used model in epidemics assumes the division 
of any fixed population size into three categories (i) susceptn- 
bles : those individuals who can catch hold of infection when an 
infected person comes into contact with them (ii) Infectives s 
those individuals who are capable of communicating the disease 
to others and (iii) Removed i those individuals who have been 
removed either by isolation, death or permanent immunity. 

Two types of models, deterministic and stochastic, have 
been investigated by researchers. Deterministic models which 
are expressed in terms of systems of differential equations for 
the numbers in each categories have been applied since the work 
of Hamer £ 1906 • But since these models do not account for 

the considerable variations which so often occur in epidemic 
situations, stochastic models are much more suitable for this 
kind of problem. Except for the case when the number of both 
susceptibles and infectives is quite large, in all other 
situations there may be important differences in the results 
predicted by both versions. 

Generally two types of situations have been studied 
in the theory of epidemics (i) simple epidemic (ii) General 
epidemic. Deterministic and stochastic models have been applied 
in both situations. 

Simple epidemic is one in which infection spreads by 
contact between the members of a community and no removal is 
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possible either by death, recovery or isolation. Therefore all 
susceptibles become infected. The contributions to this type 
of situation have been made by many. Among them are Bailey 
[] 1959, 1963, 1967a]], Kendall [] 1957]], Severo £ 1967a]] , 

Williams [] 1965a]] , Hill and Severo [[ 1969]]], Hethcote L 
and I Idle il Ql972[]. Bailey [[ 1950, 1957]] lias given a solution 
for the state probabilities p (t), where r denotes the number of 
susceptibles at any time t with a single infective initially, 
later in his paper [] 1963]] he obtained an expression for the 
generating function to' calculate the state probabilities in a 
simple way* Haskey £ 1954 [j obtained a general expression for 
mean number of susceptibles which was derived by Bailey )” 1963 * 
Williams [11965]] considered the mean and variance of susceptibles 
with one or more infectives and obtained the results for large 
populations. The case with small populations and one or more 
initial infectives has been considered by Hill and Severo [[ 1969 ]] . 
Hethcote [[ 1973]] studied the asymptotic behavior of a simple 
deterministic model in which infectives became susceptibles as 
soon as they are removed. 

The general epidemic admits removal of infectives 
from circulation, either by isolation, death or immunity. In 
any case, removals may not infect other susceptibles. Though 
McKendrick [[ 1926]] and Bartlett [] 1949 ]] had already studied 
this situation, a time-dependent solution was given by G-ani 
[[ 1965]] and Siskind [[ 1965]]. There are some new results about 
the total size of epidemic by Daniels [[ 19 67]], Downton [[ 1967a]] 
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and. Siskind Q 

an infective after recovery again behaves as a susceptible and 
that this recovery may take place a finite number of steps. 

Severo |” 1969 , 1969 b considers generalizations of some 
slochastic models and studies the case allowing for removal 
of susceptibles as well as of infect ives. 

1.5.2 DESCRIPTION 01? PROBLEMS MB RESULTS 

But in all these studies, the problem of single 

infection has been considered in a time interval (t, t + At). 

Even in birth and death processes, most of the studies are 

concerned with single birth and single death in the time interval 

(t, t + At) i.e. it is assumed that the probability of more than one 

2 

infection in this time-interval is of the order (At) . We 
consider the case when multiple infections can occur even when 
a single infective comes into contact with a number of susceptibles 
at a time. The corresponding case of multiple birth in a linear 
birth-death process has already been discussed by Doubleday jjL9 73] . 
In chapter Y we try to get insight into multiple- infect ion 
situation by considering typical cases with small populations. 

We consider in chapter Y a general case in which m 
infections can occur at a time when a infectives come into 
contact with n susceptibles initially. Por simple epidemic we 
draw the following conclusions 

Eor m < n/2 when n is even or m < n-1/2 when n is odd, 
we get a system of differential-difference equations 


1965^. Billard Q 1965^] considered the case when 
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in which the first m and the last m equations are 

written separately. But the middle n+1 - 2m equations 
can he written as one equation with r varying from 
m to n-m. Otherwise, we have to write all n+1 
equations separately. Solving the system for n = 3, 
m = 2 and n = 6, m = 2 we observe that multiple 
infection causes the process of infection to he faster 
and it takes little time for all the susceptihles to 
he infected. 'The initial probabilities tend to zero 
while later probabilities tend to 1 more quickly. 
Epidemic curves, mean and variance are also compared, 
for general epidemics, we observe that 
We write a general differential-difference equation 
in r and s, the number of susceptihles and infect ives 
respectively. Taking into account the restrictions on 

r and s and considering the probability with 

. . . . n(n+l) 

negative suffix as zero, all (n+l)(n+a+l) - — -p — — 

non-zero probabilities can be obtained. As an example 

a complete system of differential-difference equation 

is presented here for n = 3, a = 2 and m = 2 and is 

solved numerically by using Gill’s method of integration. 

Yfe also plotted the following graphs j 

(i) All probabilities of extinction for susceptihles. 

(ii) P r °bab ility that no susceptible gets 
infected and no infective is removed. 

(iii) p oq (t), the probability that all susceptihles get 
infected and all infect ives are removed* 



CHAPTER II 


OP T HAI BULOITAIION OP CPHTRO-LIB -D PIS HER IBS 

2,1 HTTROPUCTIOIT Since exploitation of common property 
fisheries had led to serious consequences even including the 
complete depletion of the resources [^Gordon ^ 1953 9 1954^ 2] > 
great attention has been paid, in the economic theory of 
exp lo itat ion of f is her ie s , to s erne s ort of c ontr ol s o t hat the 
fishery resource can be prevented from vanishing. This control 
has been applied mainly through the restrictions on effort or 
quota along with the maximization of some objectives which may 
be economic, social or both. 

Maximization of present- value of total profits over 
time has been sought by many including Beddington, Y/atts and 
fright £ 1975]], faugh and Calvo ] 19 74 ] , Sancho and Mitchell 
Q 1975, 77], Silvert and Smith £ 1977] and Clark £ 1976] . 

This chapter is divided into three more sections 
namely Sections 2.2, 2.3 and 2.4 in which we study 

(i) Sancho and Mitchell Model 

(ii) Multispecies fisheries model 

and (iii) Mathematically non- tractable models in fisheries for which 
only limiting control solutions can be obtained. 

Poliowing Sancho and Mitchell £ 1975] and a few others our 
object is also to maximize the present- value of total profits 


over time subject to constraint on fishing effort or catch, by 
the use of dynamic programming or maximum principle or calculus 
of variations in all the tnree sections. The necessity for each 
of this type of study is pointed out in the following paragraphs. 

Section 2.2 Sancho and Mitchell Q 1975^] used dynamic 
programming to obtain the optimal effort each year so as to 
maximize the present- value of total profits over time with the 
restriction on the maximum effort that can be applied. Sancho 

1975^} observe that n The results are obtained 
in terms of R, the fish population at time t = 0. If the 
maximum value of total profits is required each year, then the 
estimate of R would have to be made each year. Heuristic ally, 
this can be done by talcing the ratio of catch to effort each year 
and comparing one with the next.” 

Y/e show that the heuristic approach is not necessary 
when we use maximum principle and calculus of variations. We 
get optimal effort and optimal population as functions of time. 

As Ponte cor vo |~ 1967^] states : u Regulations specifying 
the size of capture are far more important than regulations 
determining the optimum level of effort." Y7e impose the 
constraint on fishing party not directly in terms of effort it 
can make, but in terms of duota it can harvest per unit time. 

We also find out the conditions under which our optimal 
solutions exist. This is not possible in the method adopted by 
Sancho and Mitchell £ 1975 ^ . 


and Mitchell £ 
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S ection 2 , 5 .In, this section we include the aspect that 
the area of pond or a river or a sea may contain a number of 
different types of fish commanding different prices and all of 
them may be caught. 

Section 2 . 4. In this section it is shown that limiting 
optimal solutions may be possible for an ecological system which 
may not, otherwise, be mathematically tractable. 

2.2 SMCHO AND MITCHELL MODEL 

Sane ho and Mitchell Ql975 2j have used dynamic 
programming to find optimal effort each year to optimize the 
present value of the total profits over all time. In this 
section we use maximum principle and calculus of variations to 
obtain additional results, 

2.2,1 THE MODEL . The standard model used by Sancho and 
Mitchell 1975^] , for fisheries exploitation is 

§| = aX - G(t) (1) 

where X(t) = fish population density or biomass at time t. 

a = natural rate of growth of fish 
and C(t) = catch of fish taken per unit time at time t. 

Bell and Carlson |” 19 70 took C(t) as 

C(t) = a + ]3X(t) + yB(t), (2) 

where a, (3 , Y are parametric constants and E(t) is the fishing 
effort per unit time at time t. Using (1) and (2), we get 
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= (a — p )X — a - YB(t) (3) 

The present value P of the profits is given by 

oo 

P = / e" xt QpC(t) - q.(E)3 dt, (4) 

o 

where * >0, represents the given instantaneous discount rate, p 
is the fixed price per unit of the fish harvested and q(B) is 
the cost of the fishing effort. In Clarks, Edwards and Priedlander 
[3 1973 ~[ , q(E) was taken as a linear function of E, while Sancho 
and Mitchell Q 19753 , following Holt, Modigliani, Muth and 
Simon Q I960 3 s 'took it as a quadratic function of E so that 
their mathematical problem was to maximize 

00 1 + p 

P = / e~ Ai: £p(cc + PX + TB) - (bE +kE+A)3 dt, (5) 

o 

subject to (3) and 

E(t) < Q, (6) 

where Q is the constraint on fishing effort in order to take into 
account the specified quota that may be allotted to the fishing 
party or the capacity of the market to buy fish. Using dynamic 
programming, they solve the problem obtaining the maximum value 
of P as 

f(R) = jffi 2 + BE + K, (7) 

where X(0) = R is the initial population. 

In the special case when k = 0 and z = 0, they get 
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A = + JiDi 3 . (8) 

P p + (X - 2a + 2j3)(p - pp) 

2X - 3a + 3(3 

and K = ~ Qpa + ^ (p - B) 2 Y 2 - ap ^ (10) 

They also obtain, the optimal effort 3 in this case as 

E = JK Up - - B 3 (11) 

2.2*2 PIS PUSS 101 03 EARLIER RESULT S 

Sancho and lvlitchel £ 1975 found their results in 
terms of initial population R. Talcing the ratio of catch to 
to effort each year and comparying each year with the next year’s, 
the estimate of R would have to be made each year, if the optimal 
policy is to be always followed. We show that the heuristic 
approach is not necessary since we can find for the optimal 
policy, both X(t) and E(t) explicitly as functions of time. 

Moreover the constraint that is imposed on a fishing 
party is not directly in terms of the effort it can make, but in 
terms of the quota it can harvest per unit time so that we 
replace (6) by 

C(t) = a + px(t) + yE( t ) < Q. (12) 


Initially this would imply 
a + j3R + y E <_Q. 


( 13 ) 
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If (11) does not satisfy this condition, then formula (7) 
also fails to hold and there is need for an alternate formula. 

The dynamic programming formulation also proceeds on the 
assumption that the optimal solution exists and as such the 
conditions under which (7) and (11) give the solution are not 
brought out in Sancho and Mitchell £ 1975 • We discuss the 
conditions in the next section. 


2.2.3 USE OF MAXIMUM PR IHCIPIE 

To maximize P given by (5) subject to (3) and (12), 
we form the Hamiltonian |~Pontreagin Ql96o[]» Clark £ 19 7633 

H( t,XjE,^ ) = e“ Xt { p(o+pX+YE) - ( bE 2 +kE+ SL ) } 

+ *(t) {(su-p)X-oc - YE} , (14) 

w here ip(t) is an auxiliary function satisfying 


3H _ d£ 
6X dt 


(15) 


We then maximize H as a function of E to get 
e”" X ^ (pY - 2bE - k) - Y'l'(t) = 0 


(16) 


We shall take into account the constraint (12) later. 
(Sec. 2.2,7). Prom equation (15), we get 

e PP + (a - P) ^(t) = - (17) 


Integrating (17), we get 

SP ,-xt 


V (t) = 


(X -a+p) 


e 




(18) 
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where A^ is an arbitrary constant. Substituting in (16) we get 


where 


E(t) = K 1 - Kg e ) i; J 

(py - k) (x-a+p) - ypp 


2b (x ~a+p) 


and 


A 1 Y 
k 2 = 2b“ 


(19) 

( 20 ) 

( 21 ) 


Substituting in (3) and integrating, we get 

I(t) = + .<*-*»* + Ag e'CMO* 

(X-2a+2p) * 


(82) 


where Ag is an arbitrary constant. Substituting in (5) from 
(22) and (19) we can get the maximum value of P. 


2 .2 .4 D 35 0USSI0H OP SOLUTION 

We find from (3) that if a < p , then is always 
negative and X(t) would go on decreasing and the fishery would 
soon be exhausted. 


Moreover if a < p , then x — a + p >0 and the second term 
on the R.H.S. of (22) would tend to infinity as t increases. If 
the population becomes large, our simple model 

f§ = aX (23) 


for the growth of the fish in the absence of exploitation would 
not be valid. Moreover when we substitute (22) and (19) in (5), 
we get a term 
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CO , x 

/ e (-a+p)t atj 
o 

which tends to infinity if a < {3. To avoid this we will have 

to choose Kg = 0 or A^ = 0, so that E(t) would he constant and 

a + Ylh 

X(t) would approach ~ ~ which can be negative for some values 
of parameters . 

To avoid all these complications, we shall assume 

a > p (24) 

In order to prevent the population from going to infinity, 
we shall further take 


Ag = 0 and X - a + p < 0 (25) 

Using the initial condition X(0) = R, (22) then gives 


cc + Y K- 


YK, 


R 


a-fc 


(X _ 2a + 2(3) 


( 26 ) 


s o that 

X(t) = 


a + YK 1 a + YK 1 (\ 0.0X4. 

+ ( R ~ e “ ^ 


(27) 


and (19) then gives 

E(t) = - 


AfAUU (R _ g * Tg i ) e (x-a+p)t 


"aUpT 


As t 


a + YK 1 

X(t) , B(t) 


(28) 


(29) 


If R > 


We shall assume py > k, then from (2o) and (25), > 0. 

a + YK 

TZTT 


1 


the f ish p opulation decreases to its equilibrium 
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value if R < > the fish population increases to its 

e Q. u il ibr ium value . 


a +YK. 


1 or evaluating (5), we have to find the integrals 
/ e“ Xt X(t) dt, / e~ Xt E(t) dt, / e~ Xt E 2 (t) dt 


The last integral would be finite only if 


X _ 2a + 2p < o 


(30) 


which is true from (25) . 


Yfe find from (27) and (28) that the behaviour of E(t) and 

a + YL 


X(t) is the same as t “ i.e, if E > 


T&TT 


both X(t) and 
(a + YIL) 

TO / ~L 


E(t) decrease to the equilibrium values and if R < 

then both X(t) and E(t) increase to the equilibrium values, 
This has been shown in figures 2.1, 2.2, 2.3 and 2.4. 

2.2.5 COMPARISON WITH EARLIER RESULTS 


Substituting from (27) and (28) in (5), integrating, 
using (25) and (3o) and simplifying, we easily find the maximum 
value of P to be given by the quadratic expression (7) with the 
difference that when k = 0 and fi, = 0, B has to be replaced by 
pP _ (X - 2a + 2p ) (p + 


B = 


(a-f!) 


(31) 


This difference is seen to be due to a calculation error 
in Sancho and Mitchell Q 1975]]] where equation (20) should be 
AY 2 ( p _ (3) 

*B = pP + (a— p)B - ZM 


( 32 ) 
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The second term on the R.H.S. has a negative sign and 
not the positive sign as given there. Also from (28) 

X - 2a + 26 a + yK 1 

e(o) = k x — (R - -grp--^) (S3) 

which is the same as (11) when k = 0 and modified value in (31) 
is used. 

The expressions for the maximum value of P and the effort 
required initially are the same as those of Sancho and Mitchell 
£ 1975^ when minor error in B is corrected. However we have 
also obtained E(t) and X(t) and we have not to use heuristic 
methods . 

We have also obtained the conditions 

a > p, A < su-p (34) 

under Y/hich our results are valid. 

Prom (27) and (28) we obtain 

X — 2a + 2(3 a + YL 

E(t) = K 1 - — — — (X(t) - ~y-^) (35) 

(33) is a special case of (35) when t = 0, i.e. the effort to be 
used at any time depends on the population at that time in the 
same form as the initial effort depends on the initial population. 

This is also expected. Figure 2.5 also shows that if 
a + yK 1 v 

R ^ ~ -p — = , then E < The slope of the straight line is 

positive since x - 2a + 2j3 <0. 
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X(t) - 


Fig. 2.5 
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It can also be varified that the maximum value of 


” p 

J e“ AU { p (a + px +yE) - bE } dt 


{3b 


is AX 2 (I) + B'l(I) + K, as expected. 


(37) 


8.2.6 USE OF CALCULUS OF VARIATIONS 


Substituting from (3) for E(t) in (5), we get 


E = / e - *^ |_pa + PP^ + 


(py - k) 


{ ( sw-p )X - a _ X > 


- K- { (a~p)X - a - X} 2 - H dt 
Y^ 


(38) 


Using Euler-Lagrarg* s Equations of calculus of variations 
and simplifying, we get 


wnere 


X - XX + (a-P) (X-a+p) X = 1 1# 

. , _ (PY - k)Y PpY 2 

= (x-a+p ) L + a J - ' 


(39) 


(40) 


Solving (39), we get 


X(t) = B 1 e^ x " a+ P) t + B p 


where B and B 9 are arbitrary constants. 

Now from (20) and (40 ) 

^ (py - k)Y PpY 2 

xr~^Fj = a + —^5-”””“ - 2b 


(a -P ) (x -a+p ) 


(41) 


(p Y - k) (x-a+p ) - ppY 

a + Y — - — - — —— — — — = a + yK 1 , 

2b(X~a+P) x 


( 48 ) 
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so that 


X(t) = B 1 e (»'- a +P) t + Bp 


(43) 


Since a —(3 >0, and we don’t want X(t) to increase 
indefinitely as t -*• and X(0) = li, (43) gives 


X(t) 


a + YK 1 a + YK 


alp 




(44) 


which is same as (27). 

Substituting in (3) we get (28). Thus calculus of 
variations gives the same results as the maximum principle gives* 


2 , 2 . 7 THE EIRE Cl OX C ON TROL JED QU OTAS 


If there is a quota Q fixed for the party or if the 
market demand < Q, then from (12), (27) and (28), we require 


a + 



(A - 


a + YK 


a- (3 


1 ) e (x-af£)t 


> 


(X~2a + 2p) 

+ y { K-, - — — (R 

x Y 


a + YK. 


e (x_ a+ p)tj < Q 


which gives 

(a + + (R — — (2a— X-P) e ^ a+ ^ < Q 

(45) 

If (45) is satisfied for all values of t i.e. if the 

catch finally and initially given by 

_ a + YK 1 

Y = Srp - (“ + TK i> and Z = aR + (H g-p-^Xa-x-P) 

(46) 
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is < Q, then the solutions (2?) and (28) give optimal functions* 

a +yL 

If R > , then Z < Q ==> Y < Q. 

a + YK 1 

If r < - . ^ . - 1 , then Y < Q ==> Z < Q. 

Thus the conditions for solutions (27) and (28) to he 
optimal are 

a + YK 1 

( i) R > j 2 < Q (47) 

a + YK 1 

or (ii) R < — -p— , Y < Q (48) 

a + YK 

or ( iii) R = — — , Y = Z < Q. (49) 

If these conditions are not satisfied, then the above 
solution is not optimal. If Z > Q, then the initial catch 
C(0) required by this policy would be more that what is permitted. 
We will then have to start by talcing 0(0) = Q. So long as the 
catch C(t) = Q, we have 

= aX - Q (50) 

dX 

(i) If aR > Q, then initially ^ > 0 and the population 
continues to increase and 

a + YK 1 

aX + (X - (a ~ X - p) (51) 

always remains greater than Q and there is no solution to the 

♦ 

opt imal p olicy obtained above. In this case after sane time 


tt + PX + YE = Q 


(52) 
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cannot hold because this would require E to be negative. The 
model has to be modified. 

(ii) If aR = Q and we fish at a constant rate Q, then X(t) 
and E(t) remain constant. 

(iii) If aR < Q, the population goes on decreasing till it 
becomes zero at time t^, given by 

at, 

Q = (Q - aR) e 1 (53) 

2 . 2.8 VARIATION OE E(o) AM) OE EQUILIBRIUM VALUES OE E AM) X 
WITH ‘THE PARAilEIER 'S ' ~ ' " ' ' 

(a) E(o) is given by 

p Y k(x~a+p) "d 

E(0) = ^-iUMr + - (ea - sp - x) 

ppy y(2a - 2p - X) 2ba 

2b(a~p) 2b(a~j3 ) + Y 2 

Since a > p , 2a ~ 2p - X > 0, we deduce that 

(i) E(o) increases with R i.e. initial optimal effort 
increases as the initial population increases. 

(ii) E(o) ■ increases with k if k > 0. 

(iii) E(o) decreases with p i.e. the initial optimal effort 

decreases with the increase in the selling price of the 
f ish. 

(iv) E(o) increases with b. 

(v) E(o) decreases with a. 

(vi) E(o) is independent of sl . 
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(b) The limiting value E^ of E(t) as t -*■ 00 is given by 

py - k pyp 

1 ~ 2b + 2b(a — x — p ) 

Since a-p-X > o, we deduce that 

(i) E decreases with b. 

( ii) E^ increases with p. 

(iii) E^ decreases with k if k > 0 . 

(iv) E^ increases with y . 

(v) E^ decreases with a. 

(vi) E^ increases with (3. 

(vii) increases with X . 

(c) The limiting population as t ->■ 00 j_s given by 

® -i nV^R 

x ~ = ~TTf- = ^C°: + lihr-2b + 2b^v|rxj3 

Since a-j3 >0, su-p-X > 0 , we have 

(i) X^ increases with a. 

(ii) X^ increases with p. 

(iii) X^ decreases with k if k > 0. 

(iv) X^ increases with Y . 

(y) X^ decreases v/ith a. 

(vi) X^ increases with j3 . 

2.2.9 EXP10ITATIQH OVER A FINITE H ORIZON O P TIME 

Suppose we want to exploit the fishery over a finite 
horizon of time and we want to stop just when the fishery ceases 
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to be profitable i.e. we shall stop when 

p[>+pX(I) +YE(!E)3 - ^bE 2 (l) + ItE(f) + iQ = 0 (54) 


Using equations (19) to (22) 


P { a + p 


a + YK. 


Y S A 




+ p e (x-a +P )T +B e (a_ P )T 


2b(X_2a+2p) 


+ TKl _ ^-eO-^)*) 

- b{K 1 -t^e( i -^) T } 2 -k{K 


A f e (X-Swj3)I } 

1 ~ 5TT^ e 1 


- i = 0 


(55) 


Also since X(0) = R, we get 

2 


a + YK 


£U-*fc 


i + 


Y A, 


2b(x ~2a+2(3 ) 


*+* Ag — R 


(56) 


Prom (55) and (56) we oan determine the two arbitrary 
constants A^ and Ag. In fact we may get two sets of values of 
A^ and Ag. Equations (19) and (22) then determine E(t) and X(t). 

We have also to ensure that 

p {a + pX(O) + YE(0)> - { b E 2 (0) + k E(0) + O >0, (57) 


otherwise, we do not exploit the fishery and allow the fish 
population to grow to a larger value of R. 

We may note that if we have a finite horizon of time, then 
the restrictions (24)and(25) are not necessary. 
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Even if we start with the optimal exploitation over an 
infinite horizon of time, and use (27) and (28) for X(t) and 
E(t) it would he profitable to stop when (54) is satisfied 
since any further exploitation would lead to negative profits. 

2 . 2 . 10 G OXiGLUDfflG REMARK S 

(i) If instead of (3) we use 

H = aX - dX 2 - pX - a - YE (58) 

then equation (16) of Sancho and Mitchell £ 1975]]] is modified to 
Xf = p(a+pR) + (p - | if) 2 Y 2 + {(a-P)R - dR 2 - a) || 

(59) 

If d = 0, it admits of a simple solution 
f(R) = M 2 + BR + K 

but if d ^ 0, such a simple solution is not possible. If we 
try a third degree polynomial, we get five equations to determine 
four constants and similarly for nth degree polynomial, we get 
one more equations than the number of constants. The system is 
over- determined and a unique solution does not exist. 

(ii) When we use q(E) to be cubic function of E, we get 
mathematically difficult problems by all three methods, Kapur 
and Saleem £ 1977 3* 

(iii) Yfhen we use (58) with maximum principle or 
calculus of variations we get non— linear differential equation 
for X(t) which we hope to be able to solve numerically, Kapur 
and Saleem [119773* 
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(iv) From dynamic programming solution, we can also find 
Tooth X(t) and E(t) using (11), (8) and (31) (which gives the 
correct form for (9)), we get (33) as 

(X - 2a + 2|3 ) a +-yK 1 

e(o) = k 1 - ( R - ■~icrp~”) (60) 

Since E(t) is the same function of X(t) as E(o) is of 
X(0), we get 

(A - 2a + 2p) a + yK 

E(t) = Ki _ _ —— ( Z (t) - -£^0 (61) 

Substituting in (3) and integrating the resulting 
differential equation and using X(o) = R, we get the expression 
for X(t) given by (27), substituting this in (35) we get the 
expression for E(t) given by (28). 

2.3 MUhTISPECIES FISHERIES. 

The optimality principle of dynamic programming and the 
maximum principle have been used to maximize the present- value 
of profits over all time in a fishery exploitation when the 
fishery has a number of different but non-interacting species. 

Sancho and Mitchell |“ 1975]] maximized the present- value 
of total prof its over time given by equation Q 5, Sect ion 2.2] , 
Subject to £ 3, Section 2.2] when there was only one species 
involved. 

However the area of the pond or a river or a sea may 
contain a number of different types of fish commonding different 
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prices and. all of them may be cought. Me assume that when an 
effort E is made, the catch of the ith species ( i=l,2, 3,4, . . . ,n) 
is given by 

C i(t) = a i + ^ii x l + Pi2 X 2 + ••• + Pin X n + Y i E t 1 ) 

and the selling price of the ith species is so that present 
value P of the profits is given by 


P = / e" xt 2 p.(a, + 2 B. X, + y E)-bE 2 _lcE - i 1 dt 

o Li = l 1 1 0 1 J 

( 2 ) 


v/e have to maximize P subject to 


“i 


n 


dlJ — a i X i “ ( cc ^ + 2 Pj[j E ) (i— 1 , 2 , 3 , »..,n) 

D ’ 4 

( 3 ) 

2,3.1 USE OP QPfBiALPIT PR ETQIPLE 

Let f(R 1 ,Rg, . «.,R n ) be the maximum value of P, where 
X^O) = R., (4) 


so that R^, Rg, R 3 ,..., R n are the initial populations of n 
species respectively. 


The optimality principle then gives 


f (RljRg » 



•AA 


n 


+ e f {R 1 + ACa^R^-oc^ - 2 “ y i b )j 

j— 1 

,H n + !(*-«„ - S P nj fl 3 - r n E) } 3 

( 5 ) 


9 



42 


where E Is the Initial effort and we have made use of (3) to 
find the increments in populations. 


Using Taylor’s theorem, expanding in powers of a , 
cancelling same terms in both sides of (5) and taking the limit 
as A -*■ o, we get 



This gives the following two equations to determine E and 
f (E^,Rg i • • • ? Rq) 


n 


2 (p, Y, - Y, - 2bE _ k = 0 


i=l 


i i “ i 3R 


( 8 ) 


X 

n 


n 


n 




i=l 


( 9 ) 
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The equation (8) comes because we also maximize 
f(R^, ...,R ) as a function of E. Equation (9) is obtained after 
substituting the value of E from (8) into (7). 


2.3.2 SOLUTION OB T HE EQUATIONS 


Using the quadratic nature of the problem, we assume a 
solution of (9) of the form 

n n n 


f (R^,Rg , » * . ,R^) _ 


2 2 A, . It. E. + 2 2 B. R + 0, (lo) 

3=1 i=l 13 1 3 i=l 1 1 


where A. . = A . . . 

1 J cJ >i- 

Substituting from (10) in (9) and equating the coefficients 

of R-^, Rg , * * • ,R^| R^, Rg, . . . ,R n ^ Rj^j Rp B '2 , *** jB 'n and the 

1 2 

constant term, we get g- (n + 3n + 2) equations to determine the 
same number of constants A. .*s, B.’s and G. Thus following 

J. J J- 

Bellman [1 1967 3 » a unique solution of our problem of the form 
(10) exists. 

We can determine E from (8) 


E (0) = §5- C j ^ 1 Y i fPi “ 2 (.^ 1 A ij R j + V } " (11 ) 

Equation (10) gives the maximum value of P and equation 
(11) determines the initial effort required, 

2.3.3 FUR THEIR INTEGRATI ON O P EQUATIONS 

The optimal effort at time t depends on the population 
at time t in the same way as the optimal effort initially depends 
on the initial population so that 



( 12 ) 


E(t) 


n 


i_ 

2b 1 • i vjr i 

i=l 


n 


{ Z Y,(p, - 2 Z A,, X,(t) - 2B.) - k> 


3=1 


ID 3 


Substituting from (12) in (3), we get 


ax. 

1 

aF 


n 


Y i n 


a i Z i~ a f f ^ij Z 3 ~ 2b A Y i p i + 

J— 


Y . n n 

A s s A rx(t) 

b i=l 3=1 13 1 3 


y. n Icy. 

+ ^ i=l A Yl + 2^ 


(13) 


Y/e get a set of n simultaneous differential equations 
which can be easily solved subject to 


X 1 (0) = R p Xg(0 ) = R 2 ,...,X a (°) = R n ( 14 ) 

to get X 1 (t) J X 2 (t) s . ..^(t). Substituting in (12) we get 
B(t) as a function of t. Substituting these functions in (2) 
we can get back f (RpRg, . . .,R n ) as the maximum value of P. 

The solutions of the equations are of the form 


n x-i t 

X, (t) = Z 0. . e 3 + k, 

1 j=l 13 

n Xit 

E(t) = Z f , e 3 ' + f, 
3=1 3 


(15) 


where Xp X 2 ? 1 3 j * * * s 


X n are the roots of the equation 


D - X I j = 0, 


(16) 


where D is an n x n matrix whose 


i-jth element is given by 
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For X(t) and E(t) to be finite as t we would require 

that all the roots of this equation should have negative real 
parts. The condition for this can he written by using Routtu 
liurwitz criterion or Lienard-Ohipart test vide Korn and 


Korn Q 1968“ 


2.3.4 U SE OF MAXIMUM PR BIG IP IE 


From (2) and (3) we form the Hamiltonian 


-At 


n 


n 


Ii = e { p i (a i + Z p X. + Y. E) - bE - KE -£} 


i=l 


3=1 


13 3 


n n 

+ 2 'J'.(t) { a -X-— (a- + Z p.X + Y.E)} 

i=l 1 1 1 1 i=l 13 3 1 

(18) 

The auxiliary equations are 


5H 

ax ± 


dij; . 

(i = 1,2,3, »..,n) 


(19) 


These give the n linear differential equations 


.At 


n 


n 


h p^ji + h noo - f>. u 


0=1 


= 1 31 0 


It * 
( 20 ) 


of which the solutions are given by 


’l'q('l) — ^i-2 e ^ + K i e ^ » (i— 1,2,3, . . .,n) 

w her e pi ^ * j ^n ® roots of 
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u +a 1 -p i:L -P 12 “P 13 * * * “P In 

~^21 ]i + a 2“P22 “^23 *•* ~^2n 

0 * * • • » » » • *«• 

* * • * * • ••• « • • 

“Pnl ~^n2 “Pn3 *** y+a n^nn 

The total number of arbitrary constants in (21) would 
be n. How H has to be maximized as a function of E. This gives 
, t n n 

e ( 2 p. y .-2bE-k) - 2 Y. * . (t) = 0, (22) 

i=l 1 1 1=1 1 1 

so that 

1 n , , n n u.t n 

E = jjr( S p ± Y k-e ^2 2 Y A., e 3 - 2 k. Y ) 

1=1 1 1 1=1 o=l 1 10 i=l 1 1 

(23) 

Substituting in (3) we get n linear differential equations 
to solve for ]£^(t), S g (t) , . . . ,X £1 ( t) . The solution would Gontain 
in all 2n arbitrary constants, for determinations of which we 
have t he c ond it i on s 

X^O) = R i (24) 

and that X^(t) should remain finite as t -* 

V/e solve the equations in the special case when 

hi = h’ P3.3 = 0, when j t i C^5) 

i.e. when the catch of the ith species per unit time depends on 
only its own population and the effort E. In this case (2o) 



becomes 
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atT 


•( a i - M 4 > ± - Pi P i e 


-At 


(26) 


giving 


* i ( t ) 


p , 6 . e' 
^ l r i 


• At 


a 


I ^ T - 


x + e 




(27) 


Substituting in (22) 
n 

2bE = 2 p, Y. _ k 

i=l 1 1 


At 


n 

2 

i=l 


~(a -P )t 

Y. A. e 1 1 

l l 


, ? P i t Y i 
dl ^ F x " 


(28) 


Substituting in (3) 


cGL. 


n 


tfT = ( a i“^i) X i“ a i “ 5E { ^i Y i + 


- k 


n 

2 Y A e 
1=1 1 1 


(A— a^+p j^)t 


(29) 


Integrating 


X ^ t ) = { a n . + 


a 

i T 


; Pi'ii-*) k n, i 


+ 4 

+ 2 F 


n e (X_a 3 +P 3 )t 

j=l 3 3 X- a j+ p.-a 1+ p 1 


(a-p )t 

+ M i e 11 » 


(30) 


where A.’s and M-'s are arbitrary constants. 

3 1 

Since we are assuming that > P i and * < ( a;L - f^) 
and we want X i (t) to remain finite as t -* », we get = 0 . 
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Also 

ri i={“i + 5T; J, ~s~rf~rr 


k ’i, A 

-5F- } a-- ri 


Y . n y A A^ 

+ 2 “*"”~ 

3=1 X- a. + j3. - a i + ^ 


(31) 


Those are n equations to determine the n constants 
A 1 ,Ag, . . .jA^. It is obvious that each Aj_ is a linear function 
of R 1 , Rg, . ..,R n . We then get 


y . H 

T X 


Hu 


(X-a .+§ . )t 


X- ( t) = D. + X 2 *J — e 

^ ’ 1 ^ 3=1 x- a 3 + p. - a A + fq 


3 3 


(32) 


where 


n ? Pj n< a i - x) _ j __i 


D i = <“i + sf A ~a~r-p - - x ‘ " w a t -Tj. 


i=l l r i 

n (x-a.+p.)t 

2bE = 2b E - 2 Y i A i e 

i=l 


(33) 


The maximum profit is given by 

/“ e Xt C S Pi^i+PA + YjE) - BE 2 - KE _ * 3 dt (34) 
o i=l 

Substituting from (32) and (33) in (34) and integrating, 
we get the maximum profit as a quadratic function of k v A g , . . . , iq, 
but since A^Ag,...,^ are themselves linear functions of 
R 1 ,E 2 ,...,E n , we get the maximum profit as a auadratic function 

ox R f Rg f * » * jR^ • 
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As t -* °° 


X ± (t) - O n - + 


n p. y.(a.-X) A y. 


1 ' ** Ai '4 r ^ T ~ " 


. . q r . :a Pi ( a q ~ X) 

E ( t ) L 2 Y i r — I t— r - k I3 

t=i i a i - Pi - x 


(35) 

(36) 


For n = 1, these results tally with the results obtained 
earlier. We note that 


(i) The limiting values of X^(t) and E(t) are independent 
of R ^ j Rg j • » • 5 R^ • 

(ii) When k = 0, the ultimate effort is the sum of n 
independent components due to individual species. 
However when k > 0, the ultimate effort required is 
more than the sum of the n components due to n species. 

(iii) The limiting population of the ith species is more than 
what it would be in the absence of other species. 

( iv) The effort required initially as well as at any time 
depends on the initial populations and the parameters 
of all the species. 

(v) We find 2 i (°°) = (a^ + Y q E (°°))* 

where X.(°°) and E(°°) denote tile limits of X. (t) and E(t) 

1 ax. 1 

a s t - °°. This is expected since ultimately ^=- should vanish 
in 

ax. 

3T” = “ a i ~ Y i E ^ t )* 
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Thus ultimately the population of each species is more 
than it would be in the absence of other species and the 
ultimate effort required is more than it would be for any 
individual species. As t the catch per unit time as well 

as the growth of fish of each species would be greater than in 
the absence of other species, but the rate of growth of each 
species would be greater than the catching rate of that species 
so that its ultimate population would be greater. 

2 . 3 , 5 PAR TICULAR CASE 0? TV/0 SPECIE S 

We have already solved the problem for n species. Por 
two species, we put n = 2 in all the formulae obtained earlier 
in 2 , 3.1 with an assumption that the catch per unit time of the 
ith species depends only on its population so that the equation 
( 1 ) becomes 

C-(t) = cc i + Pi + Y i E (i = 1,2,) (37) 

Assuming 

2 2 2 

f (R 1 ,R p ) =2 2 A. . R, R. + 2 2 B R + D, (38) 

1 * 3=1 L =1 10 1 3 i=l 1 1 

where A. - = A . Substituting from (38) in (9) with all the 
3-D D 1 

changes talcing into account including 1 c = 0 and x, = 0 to malce the 

2 2 

calculations simpler and equating the coefficients of R^,Rg, 

R^Rg, R 1 ,Rg and the constant term on both sides, we get the 

following six equations 


XApi = 2 (a 1 -p 1 ) A u + § (A u r ± + & 12 


(39) 
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*A 22 = 2(a 2 -j3 2 ) Agg + | (A 12 Yg + A 22 r 2 ) 2 (40) 

2AA 12 =2(a 1 -P 1 )A 12 +2(a 2 -p 2 )A 12 + | (A^Yj+A^YgXA^ Yg+AggYg) 

(41) 

2AB 1 =p 1 P 1+ 2(a 1 -P 1 )B r 2A 11 a 1 -2A 12 ag- § (YjPj+YgPg-EB^BgYg) 

( A 11 Y 1 + A 12 Yg) (42) 

2X Bg = P2P2 + ^( a 2“^2^2 “ ^ a l “ ^ a 2 "^22 

- I < WV^VW^l^V (43) 

A I) = p 1 « 1 + Pgttg - 2a - 2agBg + (l J l T i+- r - , 2 Y g - 2E 1 Y l -!iS 2, Y 2' 


Ve write equations (39), (4o) 

and (41) in the form 

(44) 

k^ A = ( Y ^ A + Yg C) 


(45) 

kg B = (Y 1 0 + Yg B) 2 


(46) 

k 3 0 = ( Y 1 A + Yg 0) (Y 1 

G + Yg B), 

(47) 


where A = k^ = b(X _ 2a^ + 2(3^) 

B = A 22 kg = b( X _ 2a 2 + 2pg) 

0 = A lg k 3 = b( X - a 1 + p 1 - a 2 + p g ) 

From (45), (46) and (47), we find 
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Substi but ing the value of A from (48) in (45), we get 
hue following two equal ion in B and C to be solved 

T 1 k 3 ° 8 + 2Y l T 2 k lh k 3 B0 + T 2 k l k 2 Be “ = 0 (49) 

0 2 + BO + Y 2 B 2 _ kg B =0 (50) 

Eliminating C, we obtain 

Q 2 = PH, (51) 

where P = B^B-kg) _ 2Y 2 k l k 2 B ^ Y 2 k 2 B - k 3 > 

Q = Y^kgCTgkgB ~ h 3 ) - Y 1 1^ (Yg B - k g ) 

R = 2Y B Y g k* ( k 2 „ k A ) fc 


We get a quadratic equation in B to be solved. This 
gives the value of Agg. Substituting the value of B in (46) 
we obtain G i.e. A^g and again substituting the value of C in 
(45) we have A^q« Similarly, from equations (42), (43) and 
(44) we obtain B^, B g and D. Thus knowing all the constants 
we can get the maximum profit f(R^,Rg). 

2.4 LIMITING- OPTIMAL CONTROn IN FISHERIES 


In this section it is shown that limiting optimal solutions 
may be possible for an ecological system which may not, otherwise, 
be tractable. We consider here the general system. As a special 
case, the system considered by Clark j~1976]] is solved completely 
for optimal equilibrium solutions. The case for interacting 
species is also considered. 
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The problem of optimal exploitation of fisheries bristles 
with mathematical difficulties. In a recent survey Clark [^1976 2] 
states 


" Our problem therefore possesses an equilibrium solution 
satisfying the necessary conditions for the maximum principle . 
However economic interpretation of this solution does not seem 
obvious. However, it appears to be extremely difficult to 
determine the optimal approach path which must consist of some 
combinations of bang-bang controls and non-equilibrium singular 
controls. We shall not make further attempts to obtain complete 
solution to the problem. The difficulty in this problem seems 
to be associated with the n on- con'cro liability of the system given 
by 


dSL 

dT 


= rx (1 - ») - q 1 Ex 


( 1 ) 


|| = sy (1 - J) - q 2 Ey, 


in the sense that x and y cannot be controlled independently ,J . 


However in the previous section 2.3 we have used both 
dynamic programming and maximum principle to obtain the optimal 
solution for the system given by 


= a x « - («! + PpX+ T^E) 

( 2 ) 

= a 2 y - (a g + P 2 y + YgE) 

inspite of the fact that here also x and y cannot be controlled 
independently. The difficulty may tnerefore really be in tne 
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non-linear nature of the functions concerned, though this nay 
ho complicated by the difficulty of independent control. 

2.4.1 BASIC EQUA TIONS 

lor discussing in a general context, we consider the 
system given by 


£f = P 1 (3 C) - EC T;L (x) - Ii-jU) 

= I' 2 (y) ~ “ H 2(y)j 


(3) 


which includes both (1) and (2) as special cases. Here 3C(t) and 
y(t) are tne biomasses of tne two species of fish, E is the 
effort, E^(X) and Fg(y) are ^ le natural rates of growth of fisn 
and EG 1 (x) + li^x) and EGg(x) + Hg(y) are the catches of tne 
two species of fish per unit time. 

If p 1 and p 2 are the two prices of the two types of 

fish per unit biomass, x is the instantaneous discount rate and 

2 

the cost of an effort E is bE + KE + A , then the present value 
of the profit is given by 

OO , 

P = / e {p (EGgx) + Hgx)) + p g (ESg(y) + Iig(y)) 

0 

- (bE 2 + kE + a ) } dt (4) 

v/e have to choose E(t) so as to maximize (4) subject to (3). 

We form the Hamiltonian j“Pontr iag in, 1963 2 
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— Aij 

H = e {p 1 (EG 1 (x) + lijU)) + Pg(EG 2 (y) + H g (y)) 

- (bE 2 + M3 + g) } + ^ 1 ( t) {P 1 (x)-EG 1 (x) - H^x) > 

+ ^ 2 (t) {Eg(y) - EC-g(y) _ £g(y) > (5) 

Using the auxiliary equations 

ail d ^l 0 li dlp 2 , rs 

W = “ TT » Ty = “ C rr 

and maximizing II as a .function of E, we get 

e Xt { Pl EG^(x) + Pl H^(x)> + ^(t) {pj_(x) - EG^(x) - H}(x) > 

dip. 

= - dt~ < 7 > 

e Xt {p 2 EG g (y) + Pgdg(y)} + *g(' fc ) {J g(y) “ S& 2^ y ) “ ^(y)* 

dt~ 

2 

dt 

and 

e Xt {p^x) + P 2 & 2 (y) — SUB - E } ~ ^(t) G^x) - * 2 (t) G g (y) 

= 0 ( 8 ) 

Prom (3), (7) and (8), we get a system of five equations 
for determining the five functions X(t), y(t), E(t), 'P-^'k) and 

♦ g(t) . 

Ve can solve for E from (8) and substitute in (3) and (7) 
to get four differential equations of the first order to solve 
for x(t), y(t), iP 1 (t) and tg( t). These can be integrated 
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numerically if x(0), y(0), ^(0) and + 2 (o) axe known. Now 
X(0) = N-p? y(0) = where R^, Rg are the initial populations. 
’J'l(O) and. 'J'gCO) axe no ^ known. Instead we have the conditions 
that x(t), y( t), 'l'^(t), 1 J'g(t) and E(t) should tend to finite 
limits as t tends to infinity. Alternatively if we want to 
e;cploit over a finite horizon of 1 time, we can use the conditions 
2(2) = 0, y(T) = 0. 

2 .4.2 Pi ffi2_I CULid ] L_ _GASE 

Equations (3) reduce to equation (2) if 
= a^iR. , 1*2 (y) = a 2 y 
Cp(x) = Y 1 5 & 2^ y ^ = Y 2 

( x ) = a-j+P-jXjflgCy) = a 2 + p2 y 
(See in Section 2. 3, 2,3 .4) . 

2.4. 3 (MR P AiJ ICUjrUi CASES 

It is easily seen that equation (7) will he independent 
of x and y if 

gJl(x) = o , Gg(y) = 0 

H^(x) = constant, fl’(y) = constant (10) 

R^(x) = constant, Fg(y) = constant 

which gives the case discussed in 2.4. 2. It is only in this 
case that we can solve for ip-^(t) and i|ig(t) separately and tuen 
get two independent equations to solve for x(t) and y(t). 
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2.4,4 (a) OP IT 3MAL E QlJILmuiai^ S 

We investigate the possibility of 


x(-t) = x 0 , y(t) = y 0 , E(t) = E 0 , 
♦ Pt) = op”, ♦g(t) = Cg e” 


( 11 ) 


as solutions of five equations (3), (7) and (8). We get 
3?l(x 0 ) - E 0 G^xp - H 1 (x 0 ) = 0 

W - E o W - W = 0 

P 1 E o a;(x 0 ) + p 1 H^(x 0 )+0 1 [:F^( Xo )-E 0 &p Xo )-Hpx 0 ): = 0 1 X (12) 

P2 E o G 2^ y o^ +p 2 1J 2^ y o^ +0 2C I '2*- y o)' J ‘o S 2( 5r o'* -H 2*' y o^ = °2 X 


p 1 G 1 (x 0 )+PgGfg(y 0 )-2bEo-l^-OiGi(x 0 )-OgG g ( y 0 ) - 0 

In general these will determine the five constants x Q ,y o , 
E , G 1 and C g but we further require x Q > 0, y Q > 0, E Q > 0 
and both x , y Q should not be zero. 

lor the particular case of (9), we get 

a i x 0 - E 0 Y i - (“i+ihV = 0 

a 2 y o “ E o Y 2 - = 0 

Ppl+ C 1 ( a 1 -P 1 ) = op (13) 

P2?2 + 0 2^ a 2-^2^ = °2 X 

Pl Y l + p 2 Y 2 - 2tE o ' E " °1 Y 1 ' °2 Y 2 = ° 
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Solving, we get 



y o ~ a" 2 -P~2 '“ a 2 + 2b 


2 , p lV x ~ a l> p 2 Y 2 ~ V q , n 

5 X — &-i + B -i X — a 0 + 8 n 


2 D 


+ p 9 Y 9 ~ k ~ 


2 ' 2 


^1 -r 

P 1 i 3 ! T 1 

X — 3^ Hh P 


2 T ^ 2 

p 2 ^2 Y 2 

a g + P 2 




Ml 


Co = 


PoP 


2-2 


'1 T- "a ~ + "JTJ 5 "2 X- a g + 8, 


It is easily seen that x Q , y , S Q are the limits of 

x(t), y(t) and E(t) as t -*■ «> and G^, Gg are the limits of 

. Xt . Xt , 

^ e , $ g e as t -* oo. 

In the general case, eliminating and Gg from the 
last three equations of (12) we get 


Pi MV + p 2 W - 2b E 0 - k 


- & h x o> 


- s g(y 0 ) 


" j _i _ E o. hiV. Hi flKL. " 

_ X - E o< x o> + E o G i( x o> + 

~ 1 = o 

_ X- Ig(y 0 ) + B 0 S'(y 0 ) + H^(y 0 ) J 

(15) 


In general this would be a cubic in E Q . For each 
positive non-negative value of E 0 > the first two equations of 
(12) would then determine X Q and y Q * 
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If b = o, we get a quadratic in E . 

o 


If however G(x) = K ± , Gg(y) = Kg, (15) becomes 


n r *n 7 ?bE >r K l P l d l (x o> t , 

p^i+Pgkg-^Eo-k- ----jr,— : - iC 2 

)+Hi(x ) 




= 0 


( 16 ) 


which is linear equation in 2 0 and we get only one positive 
value of E • 

If we consider the special case 
*l(x) = rx (1 - g), 3? 2 (y) = sy( 1 - J), G^x) = q^x, 

kg(y) = Igy J h^(x) = 0j hg(y) = 0> b=0 s Z ~ Of 

(17) 

we get from (15) 

, p l1l E o 

p 1 <i 1 x 0 +p 2 <i g y 0 -k-a 1 x 0 rrx.' 


1 _ r + - r - + E 0 q-i 


Pglg E q 

^ y o “ -r S y - 

x - 3 + T'~ + E o a 2 


= 0, 


(18) 


which gives a quadratic in E Q . Inis is the case considered by 
Glarlc Q 19 7 6 3] • 

(b) SOLU TION 


We have to solve (18) subject to 
x 
K 


r(l - r°-) = s (! - jr> = <lg B 0 


(19) 
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Substituting for x Q , y Q from (19) in (18), we get 


*( E 0 ) = (Pi^K + p 2 q. 2 L - k) ( x+ r-q 1 E o ) (x+S-q^) 

2 2 

- E 0 (Pi K r“ + p 2 i g ") ( 1 +r~q -jE 0 ) ( X + s-Q^E Q ) 
- Pl a^B 0 K(l - -UP.) (x +S -q 2 E 0 ) _ P 2 q.gE 0 L 

4oE 

(1 _ ——2) (x + r - q. 1 E Q ) 

E OW 


<{,(^co) _ oo s <J> (o) = (p^K + P 2 qL2 L - lc ) ( x+r ) (*+s), 


(j(oo) _ _co 


Also 




up- 3 -) 

!i g 


P-, Il\ (X+r) 

—nr—— * 1 % 

Pp IX (x+s) 

- — - g~ — ll*2 


x+s X+r 

. V 'ii j 


X+s X+r 

a 2 " . 


(20) 

( 21 ) 

(22) 

(23) 


Equation (21) sliows that there is atleast one positive 

value of E . Put the value of E n should be less than ~ and - 
0 0 Q-i Q-g 

since from equation (19) x Q and y Q would always be positive 
under this condition. Prom (21), (22) and (23), ?;e can be 
sure of only one such positive value of E Q lying between 

0 and when and p -i i-ht + -22^2^ ^ ^ 

q x If 9-2 x x c x, 

or 0 and when and Pqd-jK + P 2 1 2 L > k. 

2-2 q 2 2-1 
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2.4.5 SOME SOHJTICQJ 

In general when we are not able to find the optimal 
non-equili or ium singular controls, we find the equilibrium 
optimal solutions x Q , y Q , E q . If a < x R g < y , we do not 
harvest till populations rise to x Q ,y o . If either or both 
initial populations are more than (x 0 ,y o ), we do the maximum 
harvesting permissible by the market till x Q ,y o are reached. 
After (x Q ,y 0 ) have been reached, we use the optimal effort E q 
and the populations remain stable thereafter at (x Q> y )• Inis 
approach does not give the optimal solution and in fact the 
solutions may be far from optimal. We can easily compare P for 
this case with the value of P where the initial population is 

( x o» y 0 ) * 

2.4.6 .CASE .OF INTERACTING- SPECIES 

In this case, the equation (3) are replaced by 


= Fg(x,y) - EG 1 (x) - H x (x) 
= Jg(x,y) - EG 2 (y) - Hg(y) 


(24) 


and equations (12) are modified to 

F l( x o’ y o) - B o G f x o) - *A x o) = 0 
V x o’V - E 0 a 2^o) - ^Wo) = 0 


p i E 0 G i( x cA p i H it x o )+O i 




op 
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9E , 


P2 E o G 2^o^2 I %(y o ) +G 2 jd7\ 0 ,y o - E o G 2^o)“ H 2( y o> 


Pl G l( x o^2 G 2^o)- 2bE o“ k -* C l G l(- 0 ) “ °2 G 2< y 0 ) = 0 


V 

(25) 


H- 1 imir.nt iuj G., ,0-, from, the last three equations 5 we 30 1 


Pl G t x o) + fW y o> - 2b B 0 - k 

p l E o S l( x o) + p l H l x tl 


S b x o> 


5 1,1 


* - ^V,y„ +E o S I( x o) + H i‘ x o) 


0 7 ^ 0 


G 2^o) 


p 2 E o s 2(y 0 ) + p 2 0 ) 


= 0 


1 - ( tf ) S„»y„ + E o a 2 (y o) + E 2^o'> 


0 0 


(26) 


Me can solve for ic 0 ,y 0 ,E 0 from the first two equations 
of (25) and (26). If G ± (x) = y ±1 G ? (y) = Yg. (26) is a 
linear equation in E • ¥e can solve for E q from (26) and 
substitute in the first two equations of (25) to get two 
equations for solving for (x 0 ,y Q ). 

Thus for obtaining optimal equilibrium solutions or 
the limiting values of x(t), y(t) and E(t) as t -* 
of interacting species presents no additional problem 


the case 



CHAPTER III 


OPT IMAL EX PlQITAT IOh O P HIK ES 

3.1 ETEIODUCTIOH As we have mentioned in Chapter I that 
studies related to exploitation of mines have not been generally 
concerned with the effort that we apply while extracting the 
material from the mine. This aspect is studied in this chapter 

by assuming different extraction functions where the cost function 
is a function of the effort applied per unit time. 

The studies are aimed at providing the maximum present 
value of total profits over given time if the exploitation 
period is fixed and the maximum of all the maximum present values 
of total profits when the exploitation period can be varied. 

Por t lie optimal solution, the amount of the material to be taken 
out and the effort to be applied per unit time are also plotted 

against time* 

3.2 TUB MOSE L 

The model used here is 

j^=-c(t), W 

where x(t) = mass of the material in the mine at time t. 
and C(t) = mass of the material taken out per unit time at 

time t. 

Let 4 (e) denote the cost of making effort e per unit 
ttoe. It. present value (S.Y.) of total profits is given by 
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I 

P = / e At Qp C(t) - q(e)3 dt, (2) 

where X > 0 is the given instantaneous discount rate, p is 
the fixed price per unit of material taken out and I is the 
time during which the total material in the mine is taken out. 

In Clark, Edwards and Prledlander £ 1973]] , q.(e) was 
taken to he linear, hut Sancho and Mitchell ^ 1975^ took it 
to he in the more realistic form 

4(e) = he 2 + ke + £ , (5) 

where b, k, A are parametric constants and e(t) >0. 

Phus our problem is to maximize 
T 

P = f e^ £p C(t) - (be 2 +ke+£)H dt (4) 

o 

Subject to (1). 

3 .3 PJPIBB EXPRACPI02 T FU NCT IONS CONS IDE RED 

We consider the following three extraction functions 

(a) Following Bell and Carlson £ 19703 Sancho and 
Mitchell £ 1975^] , we take C(t) as linear in amount 
x(t) and effort e(t) i.e. C(t) is given by 

C(t) = a+ px+ , /e,*e(t)>0 (1) 

where a, (3,Y are positive parametric constants. 

(b) We take C(t) as linear in e(t) and Q.u-S'dratic in x(t) . 
Phis may be more realistic in practice# Now C(t) is 
g iven by 
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C(-fc) = oc+px + Yx 2 + 6e(t), e(t) > o (2) 

where oc,j3,Y and 6 are positive parametric constants. 

(c) v.'e assume that material taken out per unit time is 
proportional to the product of effort and amount of 
the material i.e . C(t) is given by 

C(t) = px(t) e(t), (-3) 

where x(t) > o, e(t) > 0 and (3 is positive parametric 
constant. 

3 .4 OPTIM IZATION POR THE FIRS T MIN Of ION 

For the first function, (Section 3.3, (l)),we have to 
maximize i ? 

P = / e £p(a+p>x+Ye) - be -ke- jQ cL^» (1) 

o 

subject to 

^ = -a - (3x - Y e (2) 

3 .4 . 1 US E OP MAXIM U M PRI NCIPLE 

To maximize P given by (i) subject to (2), we form the 
Ham ilt on ian £p on try ag in Q 19 63 [] 2 

H(t,x,e,<J») = e Xt {p(a+px+Ye) - be 2 -ke~ - 

4>(t) (a + (3x + Y e), (®) 

where ty(t) is an auxiliary function satisfying 

9H cW -.('4) 

9 x 

We also maximize H as a function of e to get 



67 


Substituting from (8) and (12) in (1), we obtain 


^ k 2 . -AS. B 2 BT .. 

P = -- (1 - e ) + — - (e K - 1) + 


8 


(*+P) 


2 


} 


B 


'2 ,_(X+2p)l 


( X+2p ) 


{e' 


1 } 


where 

2 

Ag = pa - ppL^ + K^(pr - k) - bK^ 

Bg = pB 1 + {2b K 1 B 1 - (pY - k) B 1 > 

Og = p c 1 

b(A+2p) 2 B? 
and Bg = — - — — g— — -= 

3.4.2 U SB OP CALCULUS OP VA RIATIONS 


(16) 

(17) 

(18) 

( 19 ) 


( 20 ) 


Substituting from v (2) for e in ( 1) we have to maximize 


T 

P = / £p(-X) - -^p (cc+Px+x) 2 + ~ (a+px+x) - jQdt (21) 

o Y 

Using Euler-lagrange * s equation of calculus of variations and 
simplifying we get 


x - Ax - P (A +(3 )x + 


(A+P) (pY — 2ka) - Ap 


2b 


= 0 


( ( 22 ) 


Integrating (2a) with the boundary conditions x(0) = B 
and x(S) = 0, we get the same equation for x(t) as in (12)* 
Substituting this value for x(t) in (2)» we get the same 
expression for e(t) as in (8) and putting these two values in 
(1) and simplifying we get the same expression for P as in 
( 16 ). 
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3.4.3 D ISCUSSION OP SOLUTION 

The expressions (8) , (12) and (16) are true for any 
arbitrary set of positive values of parameters. We study their 
variation with respect to t and 1 for the following set of 
parameters 

a = 1.0, P = 1.0, x = .1, Y = 4.0, k = 1.0, p = 500.0, b = 1 
and x(o) = 100.0 

Pig. 3.1 snows the variation of x(t) with respect to t with 
T = .035, .090, .290, .350 

Pigs. 3.2, 3.3, 3.4 and 3.5 show the variation of e(t) with 
respect to t for 

I = .035, .090, .290 and .350 
respectively. 

Pig. 3.6 shows the variation of P with respect to 1. 

Pig. 3.1 shows that with the increase in I, the rate of 
decrease of x(t) slows down which is obvious. If we take I 
very small, x(t) decreases very rapidly almost linearly. 

Pigs. 3.2 and 3.3 show that for exploiting the mine 
in a short interval we go on increasing the effort while on 
the other hand exploitation in a large interval leads to 
continuous decrease in the effort Pigs. 3,4 and 3.5. 

Pig. 3.6 shows that if we exploit the mine in a time 
interval which is smaller than a certain limit we go through 
a loss. But exploitation beyond this limit leads to the 
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profit and this profit goes on increasing with I. The figure 
also illustrates that we do not get any additional profit if 
we exploit the mine in a very large interval. 

3.5 QP L IMI ZAHON FOR THE SECOND PUUCfl OI 

lor second function (Section 3.3, (2)), we have to maximize 
I 

P = / e Xt £p(ct+px+Yx 2 +6e) - he -he- 5~] dt (1) 

o *~ 

subject to 

■~2=_[^a + px(t)+Yx 2 (t) + 6e(t)H[ (2) 


3.5.1 USE_ OP CALCULUS OP VARIAT IONS 

Substituting for e(t) from (2) in ( 1) , we get 

P = f e X-b f-px - (a + px + Yx 2 + x) 2 
o 6 


+ ~ (oc+px + yx 2 + 5) - dt 


( 3 ) 


which is to be maximized. Using Euler-Lagrange ' s equation 
of calculus of variations and simplifying, we get 

x — Ax — 2Y 2 x 3 - Y(X+3p)x 2 + (k6-2ba)-P (X+P ) U x 


(X+p) (16 - 2ba) - xp6' 


“2F 


= 0 


( 4 ) 


which has to be integratipgf subject to initial and final 
conditions 

x(0) = R and x(L) = 0, 


1 


where E. is the total amount of material in the mine. 



3.5.2 PROCEDURE POR SOLVING THE PROBIBM 

We have got to solve a second order non-linear equation 
in x which cannot he analytically integrated. We have here to 
solve a two-point boundary value problem. We integrate (4) 
subject to the initial conditions 

x(0) = R, i(0) = v Q (6) 

and continue the numerical integration till the material in the 
mine is completely exhausted i.e. ,till such time T where 

X(T) = 0 (?) 

Dor each value of v , we get in general a unique value 
of I. We try various values of v Q till we get the value of T 
given to us. 

In practice we would like to treat T also as a control 
variable. We calculate x(t) and ^ for various values of v Q 
and for each v 0 > we find 33 • 

dx 

Knowing x(t), ^ and 33, we can find e(t) from (2) and 
then find P from (l) . Por each v D there is a unique value of 03 
and unique value of P . We can choose v Q so as to maximize P. 
Then we can choose the optimal period of exploiting the mine 
so that if choose this period and use the optimal mining policy, 
we shall get the maximum of all the profits. 

3.5.3 DIMEE'S IOWIESS POSH 


Making the substitution 



equation (4) reduces to two first-order simultaneous differential 


equations in x and v as follows 
dx 

= v 

XV + 2y 2 x 3 + y( x +3P)x 2 - Zl (Jc6 - 2bcc) - P(X+P)^x 
-U+P) (3s& - 2ba) - xp6 2 Co) 

— QV, \ - J 


Using the non-dimensional variables 


a 


, V = f and x = tp. 


the equation ( 9 ) reduces to non-dimensional form 
dX 


dx 

dV 

dr 


= V 


AV 


+ BX 3 + (0 1 + 3Gg) X 2 - DX + 1, 


where 


„ 2 2 
2 y a 


— r" » D 


Ya 

C 2 = “ 2 T > D = 

^ p 


YXa 

T~ » G i = TT 
Y(k6 - 2ha) - bp(X+p ) 


up 


2 " 


p * (j&_r-g££l 

2bap 


Equation ( 2 ) be c ome s 

V(t) = g = -1 - X - C 2 X 2 - E(t), 


(- 10 ) 


( 11 ) 


( 12 ) 


( 15 ) 


where E( t) = e(t) . 

3.5.4 A IHJMBftlCA L IILUSTRATIQI 

We integrate (ll) for the following values of parameters 
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A = 0.1, JB = 2.0, C 1 = 0.1, C 2 = 1.0 

(14) 

D = 2,1, 1 = 50.55, 
and for the initial value 

X(0) = 100.0 ( 15 ) 

From (13) 

V(0) = - 1 - 1(0) - X 2 (0) - E(o) 

= -10 101.0 - E(o) (l 6 ) 


We integrate 
dX 


37 


= V 


3“ = 0.1V + 2. OX 3 


+ 3. IX - 2. IX + 50.55 


subject to 

X(0) = 100, 7(0) = -10102, -10201, -l050l,-Hl0l, 

... (17) 


for which E(o) = 1, 100, 400, 1000. 

lor each of the four cases, we 


to get the following table 



E(0) 

HO) 

(a) 

1 

-10102 

Cb) 

100 

-10201 

(0 

400 

-10501 

(*) 

1000 

-lUOl 


find the optimal profit 


pi 

P 

147 

98970.60 

040 

67495 .40 

,025 

5980.64 

.017 

-98044.00 


In case (a) initial rate of exploitation of mine is 
least, the time taken to exploit the mine fully is largest and 



the present value of discounted profits is maximum. As the 
initial rate of exploitation increases, the total time of 
exploitation and the profits decrease fast and a stage is 
reached when the profit becomes zero and if we exploit the mine 
still faster, the exploitation of mine would lead to a loss, though 
due to use of calculus of variations, the losses should he 
minimum. 

Fig. 3.7 shows the variation of X(t) with r for the 
four cases. 

dX 

I l 'ig. 3.8 shows the variation of ^ with t for the 
four cases. 

Figures 3.9, 3.10 and 3.11 show the variation of E(t) 
with t in the first three cases. 

In the first case, the effort starts with initial 
value 1 and increases steadily to the value 130 over the time 
interval .147. 

In the second case, the effort starts with the larger 
initial value of 100 and increases rather guickly to 1420 over 
the time interval .040. 

In the third case, the effort starts with the largest 
initial value of 400 and increases still more quickly to 2940 
over the time interval .025. 

In this example it appears that slower exploitation of 
mine over a long period is more profitable than quicker 
exploitation over short time. 




<K < — 



fig. 3.7 Pig. 3.8 



Fig. 3.9 
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3 . 6 QPTI MIZA1‘I0IT ICR THE THI R D MCTM 

lor the third function £Sect ion .3 .3 , , (3)), our problem is 

to maximize 
!T 

1 = / e^ j_pj3 X(t)e(t) - be 2 (t) - ke(t) - Ji^dt (1) 
o 

subject to 

H = -P *(10 e(t) (2) 

3.6.1 U SE 01 CALCULUS OF VARIATIONS 

Substituting for e(t) from (2) in (1) , we get 

p ELp(-*) ~ + jjl ” O dt » ( 3 ) 

which is to be maximized. Using Eulr-Lagr ange ' s equation of 
calculus of variations and simplifying, we get 

xx ~ x 2 - xxx + x 2 - x 3 = o (4) 

which is to be integrated subject to the initial and final 
conditions 

X(0) = R, X(l) = 0, <05.) 

where B. is the total material in the mine. 

3.6.2 H OW 10 SOLVE THE PROBLEM 

We have here to solve a second order non-linear 
differential equation in x which can not be integrated analytica- 
lly. We have got here to solve a two point boundary value 
problem. We integrate (4.) subject to the initial conditions 
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*(0)=*, if| = v (6) 

0 

and continue the integration till the material in the mine is 
supposed to be exhausted i.e, till such time 1 where 

*(I) = 0. (7) 

lor each value y q , we get in general a unique value of 
1. As in earlier cases, here also we would like to treat I as 
a control variable. Me calculate x(t) and for various 
values of v Q and for each v Q , we find 1. We can choose the 
optimal period of exploiting the mine so that if choose this 
period and use the optimal mining policy, we shall get the 
maximum of all the profits. 


3.6.3 D IHENS IONIESS FORM 

Making the substitution w = v, equation (4) reduces 

to tv70 first order simultaneous differential equations in x and 
v as follows 



dv V 2 llr 2 

^ = r + U “SF * + ~£1T *• 


( 8 ) 


Using the non-dimensional variables 


V = jjf , T = tp 

equations (55) reduce to non-dimensional form 


ax _ 

dr ~ 


V 


d J _ vf. 

d i ~ HL 


+ AV _ B% + Gjfe 2 , 


(9) 


( 20 ) 
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where A = p- , B = g^p» 

XP 

and C = • 

Equation (49) gives 


( 11 ) 


7(t) = -e(x) *(t) (12) 

3.6.4 I TOI.1ER ICAL ILLUS mATION 

We integrate ( lo) for the following values of the 


parameters 

A = 0.1, B = .25, 0 = 5.0 (13) 

and for initial value 

*( 0 ) = 10.0 ( 14 ) 

Prom (59) 7(0) = -e(o) X(0) 

= -10 e(0) i( 15) 

',/e integrate 

®=V, ^=“+ ,1V- .25X+5X 2 (16) 

subject to 

*(0) = 10.0, 7(0) = -loo, -150, -200 (17) 

for which e(o) = 10, 15, 20. 


Por each of the three cases, we find the optimal profit 
to get the following table : 



e(0) 

v(o) 

pi 

P 

(a) 

10 

-100 

1.177 

933.001 

0>) 

15 

-150 

0.375 

904.343 

(c) 

20 

-200 

0.254 

879.393 
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Our condition on e(i) is to be positive for all times 

till the mine is exhausted. For all values of J Q 
‘teaoifi i 

which are ifoss than -100. 0, e(x) comes out to be negative 
some time before the mine is exhausted. That is why we start 
with V(o) = -100.0. 


In this example, in case (a) initial rate of exploitation 
of mine is least, the time taken to exploit the mine till the 
material left in the mine is 1 percent is largest and the 
present value of discounted profits is maximum. As the initial 
rate of exploitation increases, the total time of exploitation 
and profit decrease. 


Fig. 3.12 shows the variation of 3£(t) with r for the 
three cases. 

Fig. 3,13 shows the variation of ^ with t for the 
three cases. 

Fig. 3.14 shows the variation of e(x) with x in the 
three cases. 


In the first case, the effort starts with initial 
value lo and decreases steadily to the value 2.74506 in the 
time interval 1.156 then increases to the value 2.74540 till 
the material left in mine is 1. 


In the second case, the effort starts with the larger 
initial value 15 and decreases rather slowly to 11.7706 in 
the time— interval .288 then increases to the value 11.8199. 
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\ V=~200 
18 I 



T *“* 


Fig. 3.14 


i 
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In the third case, the effort starts with the largest 
value 20 and decreases still more slowly to the value 17.8161 
in the time interval .175 and then increasesto the value 17.891. 

In this example also it appears that slower exploitation 
of mine over a long period Is more profitable than the quicker 
e xp lo it at ion over sh or t _ t ime , 



CHAP TEH IV 


BIRTH AND HEATH PROCESSES fflTjQn &gA.™! 

This chapter runs over two sections 

4.1 Continuous-time linear birth-death-emigration process. 

4.2 Discrete-time linear birth-deat^immigraticn-emigration 
process . 

In continuous time birth-death-emigration process, we 
obtain the probability of extinction numerically, in terms of 
the solution of an Abel’s integral equation, which has been 
left out by Getz £ 19751] and pointed out by Kapur C 1979s O * 

In discrete-time linear birth-death- immigration 
emigration process, a discrete version of the continuous time 
birth-death-immigration-emigration process is given and solved 

numer ically . 

4 . i GCHTnrU OIJB- TIllE^ LIHEi ffi 

4.1.1 IHTR O DDC 1 M. Let n denote the size of population andtfn.t) 
the probability of betas » Persons ta the system at any time t.Iet 
the probability of a birth, a death, an tamferatiah and an 

emigration be nl it * O(it), nvat + 0(4t), vit + 0(4t) and 

/ . \ o. • lv Y fe obtain the following system of 

a At + o( At) . respectively. uuw " 

differential- difference equations 
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P (n,t) = C (n+l)y+a 3l( n +l» n(y +X )+a+ P(n,t) 

+ |I (n-l)x+ v] P(n-l,t) n=l, 2,3, . . . (1) 


P*(0,t) = (y+a)(p(A,y _ vp(o,t) 


•Defining the probability-generating function 

oo 

G(s,t) = z P(n,t) s n , (3) 

n=0 

multiplying (1) by s 11 and (2) by s° and summing up for all values 
of n, we get the partial differential equation 

§§ = ( X s-y)(s_l) || + (v - |)(s-l) G + =11 P(0,t) (4) 

Getz (3 1975,76]] obtained (4) without the last term 
on the right hand side and solved it, but the omission of this 
term containing the probability of extinction P(o,t) made all 
his calculations wrong. This error was pointed out in ICapur 
[2 1979a 3] and was accepted by Getz [3 1978 3* Ihe P ro ^ em of 
solving the complete equation (4) however remained. Por the 
steady case this was accomplished in Kapur |jL978a, 78b,78c,78d, 79dj 
Por the non- steady case, we obtain the solution below for the 
case v - o, cc < y in terms of the solution of an Abel’s integral 
equation. Another solution in terms of incomplete beta functions 
has been given elsewhere Kapur j~ 1979b, 19 79c 3) • As expected 
the agreement between the two solution is quite good. 
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4 * 1 * 2 BASIC HilEGRAL EQUAT ION 

If n Q denotes the number of persons initially in the 
system, the solution of (4) is found out in Kapur |_ 1979 a 3 • 

let P(o,t) = p(t) - (5) 

Then equation (4) becomes 

- (Xs-y)(s-l) || + (v - |)(s-l) G + p(t) (6) 

I he auxiliary equations are 

dt _ do dG / r/\ 

r- - ixs-mki-s 7 " (v _ |)(' s _i) G + all=ii P (t) 

s s 

One intermediate integral is 

=JL e( X - u ^ = iV-H z = IC (8) 

x— S JU-S 

lor the second intermediate integral 


_J.s _ _ d G _ 

rxsiyyrn'sT (v _ - )(s :T)g + 


( 9 ) 


we get the equation 


dG 

cCz 


(yv - ccX) z + K(v-cc) 
(yz + K) (X z + K) 


G 


- ^ 


( 10 ) 


Integrating and using (8), 
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a 


(l-w) 

is t 


fx-y 

Lrfc^ 

~ Z J 

|{T j 

Lt-~!j 


-vA 


a 


z a _1 

|/ f(z)(yz+K) y (XZ+K)“ V//X dz 


K = z 

1— s 


+ z) 


( 11 ) 


where ^(u) is an arbitrary function of u. If initially there 
ere n Q persons in the system, we get 

a v 


n 


~ X-y ~ 
_ 1-s 


y - r 2. 
y 

s = ip 


1-s 


(12) 


giving 


*(y) = 


n + - 
o ^ y 


(y + x) 


% +V/T 


(13) 


So that 


a v a v a 


G ( s ji) (iTs ^ W zV 


X s v = 


a 


z a __ i y 

f(z)(y 2 H-K) y (xz+K) X dz^J 


K = 


Xs r .. M, z 
1-s 


n o + u 


< a cr IL "»»> 

. _ n +v7x 

2 + x) 


( 14 ) 


Putting s = 0, we get in (14) and simplifying, we get 
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x 


/ f(z) (nz 
1 


vx) y 


£ - 1 


(Xz - px)“ V/ ^ dz 
n + a/p 


_ i iiL.- ji£) — __ 
a & +v/x 

(X - px) ° 


This is Volterra's' Integral equation of the firs 
whose solution should determine P(0,t). 

Putting v = 0, we get Abel's integral equation 


/ f(z) (z-x) 
1 


n +k 

k-1 dz _ _ 1 U^LL—. 

az> ~ F n 

(K-x) 0 


yfhere k = — and K = X/p 

v 

4 . i . 3 SOLUTION OP INTEGRAL EQUATI O N W X < p 

It is known My skis £ 19 75 3] that when 0 < a < 1 
solution of 

X 

/ _i- f ( z) dz = p(x), C = o 13 

b (x-z) a 
is given by 


*U) 


sin wa r 

~^~l 



When X < p, z >1, x > 1, l<z<x 
Clfc) 


Equation can be written as 
x 


n 0 +k 


r f(z)dz 1 (x— 1) 

{ T~X^ = E '^"*r 

1 ( x-z ) (x - K) 


(15) 

t kind 

( 16 ) 

(17) 

the 

(18) 

(19) 

( 20 ) 


(El) 
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When 0 < k < 1, its solution is 


<L ilrll 


f(x) = a a . aL V isl | / 


n o+ lc 
n~ 


* dy (y-K) 0 


(x-y) 


ir 


iy 


sin x( 1-k) 


kit 


n. 


[ ( n o +k > { 


x (y-1) 


n +k-l 
o 


n 


iy 




n 0+ k 


n 4-1 

(y-K) 0 (x-y) 


1 (y-K) 0 (x-y) 
dy 


k 


(22) 


The two integrals on the R.H.S. have to he evaluated 
nuraer ically . 

4.1.4 SOLUTI ON Ok THE INTEGRAL EQTJATICHW HM X _>_L 

In this case x < 1, z < 1, x < z < 1, K > 1 and the 

solution is given hy 


f(x) = f (n 0+ k) / 


n Q +k-l 


-O'"' ° & V 

x (K-y) (y-x) 


dy 


- n °£ 


13 


n 0 +k 


(K-y) 


n 0+ i 


(y-x) 


TV dy 


] 


(25) 


The two integrals on the R.H.S. have 


to he evaluated numerically. 


4.1.5 SOLUTION OR THE INTEGRAL EQUATM _WM X-y 

In this case we cannot use the substitution z=e 
Equation (4) can he written as 
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0G 

Tf 


»( S -1 ) 2 || + (V - |)(8-1)S + 2%U fi(t) 


(24) 


Its solution is 

O oc/p t 1 , 

G ( Sj t) - 3 —— - = / a p(I) (vl + A- - ut)- v / p 

( S— 1 ) — O s “ 

a 

(vl + - vt + l) v dl 

+ C—J - U t)- v/ '' (=lj - Mt + 1 )" (- 3— 1 + I ) 1 " 0 

S=t - ^ 

(25) 


Putting s = 0, we get 

"tj 

a / p(T)(pt+l-pf )“ v / y (vt - pi)' 1 ~ dl 
o 

= «(1 + t a/i ( r n_ f) n o (26) 

which is again a Vol terra’s integral equation of the first kind* 


When v = 0, we get 



eigM— 

(t-i)^ 



(27) 


Its solution is 


.. t * 

’-fhl r 

i 1 

0 (t-i) 


. _ n -1 

^o(rat) 

~ 


— rr d£ 

( 1+Vl) 
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s in -rc( l~k) 

Ttk 



kT' 


,k~l 


(t-T) 


T ( 


VI \ 
1+yT' 


n 

o 


dT 


, n -1 
t i lc (uT) 0 

+ n Q ix f — ■ — dT 

0 ( t~T ) k ( 1+ yi ) 0 


(28) 


The integrals on the R.H.S have to be evaluated numerically, 
4 ,1.6 ra.'LERICAl INTEGRATION 


Vfe use the method of L.V. Kantorovich, Demidovich ["] 1976]] , 
for isolating singularities to calculate the approximate value 
of the integrals as follows 

/ 4lll.lL. _ / 

J i V — J 1 V 


o (t-T) 1 ^ o y 1 ^ 


_ _ 4!ill £ 


k+1 


k+1 


+ 




k+2 


t 

E+2 


TJ k+3 +J o 


The integral on the R.H.S. is a proper integral and 
is computed by Simpson's formula. The results are illustrated 
in the following table. The numbers in brackets give the values 
obtained by using the expressions in terms of incomplete beta 
functions Kapur [] 1979b, 1979c]] and the tables of incomplete 
beta functions Pearson []1968]]. The small differences are 
due to approximate numerical integration in one case and 
linear interpolation in tables the others. 
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4 • 2 PIS CR ET E- TIME BIRTH- D EATH- iMlIGlUTIOI-EMI&RAIIOfkPiiOGE S S 

4,2.1 IH3_ BAS IC EQIJAf IQ1IS M B THEIR SOnUHOH 

Ihe basic system of differential-difference equations 
for the continuous- time B-D-E-I process is given by 


P ’ (n, t) = P(n~l, t ) [_(n-l)x+vlj - P (n,t) £ n(y+X)+a+ v^J 


+ P ( n+ 1 , t ) (^ ( n+ 1 ) y + a“| n — Ij2j0y«.» 
P'(0,t) = - vp( 0 ,t) + (y+a) P(l,t), 


( 1 ) 


where P(n f t) is the probability of there being n persons in 
the system at time t. 

We now consider the discrete version of (1). viz. 

P (n, t+1) = P(n~l,t)|” (n-l)x+vH+P(n,t)C l-ny-nx-a-v^ | 

+ P(n+l,t ) (n+l)y+a]3 ri = 1,2,3, .. . 

P(0,t+1) = (l~v) P(0,t) + (y+a) P(l,t), 

We can regard (2) as the basic equations for a 
discrete B-D-.E-I process. Given the probability distribution 
at time t, (2) will enable us to find the probability distribution 
at time t+1, so that if we know the probability distribution at 
time 0, we can find the probability distribution at future times. 

If there are n Q persons in the system at time o, 

P(n Q ,0) = 1, P(n,0) = 0 when n ^ n Q (3) 


1 — a- v ' 
X+y 


( 2 ) 



l e t n 0 = 3, then by using the above equations, we get 


2(0,0) = 0, P( IjO) = 0, P(2,0) = 0, P(3,0) = 1, p(4,0) = 0 
P(5,0) = 0, . . . 

P(0,1) = 0,P(1, l)=0,P(2,l)=3y+a,P(3,l)=l-5i-3A~a-.v, 

P(4,l) = 3X + v, P(5,l) = 0, P(6,l) = 0, ... 

P(0,2)«0,P(l,2)=(2p+a)(3y+a),P<2,2)=(3P+a)(2-5p-5A-2cu2v), 

P ( 3 , 2 )= ( 1-3P-3 X_a_ v ) 2 + ( 3y +a ) ( 2X+ v)+ ( 4 P +a ) ( 3X+ v ) , 

P (4, 2 )=(3 X+ v)( l~4y— 4X —tx— v)+ ( 1— 3y— 3X — (X— v) ( 3 x+ v) , 

P(5,2)=(3X+ v)(4X+ v)j P(n,2) = 0 if n > 5 __ 

We can similarly find P(n,t) for t = 3,4,5,... and 
varying from 0 to 3+t. 

It is easily verified that 

oo oo 

2 P(n, 1) = 1, 2 P(n,2) = 1 

n=0 n=0 

OO 

(n), i = 2 nP(n,l) = 3+3X..3y-a+v 

1 n=0 

OO 

(n), p = 2 nP(n,2)=(3y+a)(4-4x-8y~v-3a) 

+ (3X+ v)(8-8x-*16y-3 v-5a) 
+3(l-3X_3y-a-v) 2 



97 


4.2.2 SC MB RESUL TS ABOUT TH E PROCES S 


(i) If at t = o, there are n Q persons in the system, then 

P(n,t) = 0 if n > n Q + t or n < n Q - t (8) 

(ii) Thus the non-zero probabilities are ; 

P(n Q ,0) 

P(n 0 -l,l),P(n 0 ,l),P(n 0 +l s l) 

P(n 0 -2,2),P(n 0 _l,2),P(n 0 ,2),P(n 0 +l,2),P(n 0 +2,2) 


■2(0, n^) ,P( l,n Q ) ,P(2,n Q ) , • . • • ,P ( 2n^ ,n o ) 

(9) 

P(0,n),P(l,n),P(2,n) , ,P(n 0 +n,n) 


where the sum of probabilities in each row is unit. 

(iii) The number of non-zero probabilities upto time 
n (> n ) is 

= 1+5+ 5+ . . «+( 2n Q +l) *~| + ( 2-^1 q+2 )+ ( 2n^+o )+ . . . 

... + (n 0 +n+l)[] 

= g [^n 2 + 2nn Q + 3n + n Q - n^ + 2 ^ (10) 

(iv) Using (8), equation (2) gives 

P( t+l+n 0 , t+l)=P( t+n Q , t) £ ( t+n Q )X + v ]] 

= P( t+n Q -l, t~l) Q t+n 0 -l)X+vj [11«i 0 )x + v] 
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= L ( t+ %) x+v I] C ( t+n o -1 ) x+ 
... Oo X+V U p ( n o’°)’ 


( 11 ) 


so that 


n o +t~l 


P (t+n , t) = n (eh + v) 

^ *wv VXO 


( 12 ) 


m=n. 


and v . 


T/e may note that this probability is independent of cc 

(v) Multiplying (2) by n and summing for all values of 
n, we get 

(n) (1 + x _ p) n( t ) + v - a + a P(0,t), (13) 

'C-4" 1 

which can he compared with 

H = ( X -y)n+v-cc + a P(0,t), ( 14 ) 

obtained in the continuous case [__Kapur [ 19?9a] ] • 

(vi) Multiplying (2) by n and summing for all va 

of n, we get 


(n 2 ) ^ 1= ( 1+ 21 -2 v) (n 2 ) t + ( 2v ~2a+X+ v) n 

+ v + a - a P(0,' c )) 

which can be compared with 

E(n 2 )=2( X- u )E(n 2 )+(2v-2a+x+y) n 


(15) 


(16) 


+ v + a - a P(0,t) , 
obtained in the continuous case [^Kapur C 1979a 33* 
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(vii) Getz £ 1975, 19 76 obtained (14) and (16) without 
the terms involving P(0,t). ihese made the 
calculations of (n)^. and (a )^ simpler, but the 
formulae obtained were incorrect, Ve cannot also 

g 

use (13) and (15) to obtain (n)^. and (a )^ unless 
we first obtain seperately the probability of 
extinction for all t. 

4.2,3 GR APHS 

(i) In Figures 4.1 - 4.3, we have plotted P(n,t) against 
n for t = 1, 11, 21, 31 for X =0.3, y =0.1, 

a = .05, v = 0.02 and n Q = 3,4,5 (17) 

(ii) In Figures 4.4 - 4.6, we have plotted P(n,t) agaiist 
n for t = 1, 11, 21, 31 for X = 0.03, y = 0.1, 

a = 0.0, v = 0.0 and n Q = 3,4,5. (18) 

(iii) In Figures 4.7 - 4.8, we have plotted P(n,t) agaiist 
t for x = 0.03, y = 0.1, a = 0.05, v = 0.02 and 

n Q = 3,5. (IS) 

(iv) In Figures 4,9 - 4.10, we have plotted n(t) and 

2 

,0 (t) against t for 

X = 0.03, y = o.l, a = 0.05, v = 0.02 and 
n Q = 3,5. 


(80) 
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m IdULTIPia IHBBCII O h Bt EP XaSKIGS 

5. 1 IK aEODU GIIOK : As we have mentioned In Chapter I, all 
the work in epidemics whether on the theory development or 
its applications has been based on the assumption that only 
single Infection take place in the time interval (t,t + At). 

In other words, we have assumed that the probability of two 
or more infections in this interval is an infinitesimal of a 
higher order than the probability of a single Infection in 
this interval. I he case with multiple births in birth and 
death process has already been considered by Doubleday C 1973 J , 
Kapur Q I9 79e]]] and Kapur and Kumar £ 1978^ • We assume here 
in this chapter that more than one infection can take place 
in (t,t + At). Our attempt is to get insight into such 
situations. lo begin with, we consider very simple models in 
epidemics based on Ba iley C 1957;]. 

The first model is a simple epidemic model which 
describes the situation when susceptibles are getting infected 
with no removal. The second is general epidemic model which 
takes into account that infect ives are being removed either 
by death, isolation or permanent immunity. Both deterministic 
and stochastic versions are considered and compared. 
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5 ■ 2 SI MP IE BPIPEIIIC MO DEL 

5,2.1 PE HEEL IK I S 'I ' I C_ MO IEIi : Let x and y denote the number 

of susceptibles and infect ives respectively at any time in 
the total population size n+1 with x(0) = n and y(0) = 1* 
let p , 2 pg, 3Pgj.,. ? ... mp m denote the infection roles of 
there being one person infected, two persons infected,..., m 
persons infected at a time. The number of new infections in 
time dt is (p 1 +2p 9 + m6 ) xy dt. The differential equation 
which describes the process is as follows 

~ p ± xy - 2p 2 xy _ . - mp m xy (m < n) 
t w it h x + y = n +• 1. 

Changing the time scale to t = tp we get 


— — (l + 2 g 2 + SGg + ... 4- mS^) xy , 

p. 

where (i = 2,3,...,m) 

x p ^ 

Substituting the value of y, we obtain 

~~ = -( l+2e 2 +Se 3 +. . .+me ri ) x(n+l-x) 

with initial condition 
x(0) = n 

Integrating (3), we get 

n(n+l) 

X ( T ) = — — 

n + e 


( 1 ) 

( 2 ) 


( 3 ) 


( 4 ) 


( 5 ) 
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when there is only single infection i.e. e 0 = e* = 0 = 0 

iC o m 

then 


x(t) 



n + 


( 6 ) 


We observe that at any time t, the number of susceptibles 
are reduced since j3^ in the case of single infection in replaced 
by + ... + mp m in the case of multiple infections. 

cl"x 

for the epidemic curve, we find out - to know the 
rate at which new cases occur 


^ n(n+l) (l+2e 2 +3e 3 +. ..+me m )e 


( 1+26^+36^+. . »+me m )(n+l) T 


dt 


( l+2e r>+ . . »+mG^) (n+ 1 ) t ^ 


{ n + e 


This curve attains its maximum at t = 


(?) 

— — ISSJSL 

( i+2e 2 + • • 


and - 


dx ( 1+2S 2 +3 g 3 + * « .4-mG^) 


dT 


Epidemic curve given by (7) is plotted for n = 6, m = 2 

in figures 5.1, 5.2 for e = 0.2 and 0.5 where it is compared with 
the corresponding curve obtained for the stochastic case obtained 
in Section 5,2.3. 

5.2.2 SfOCi UST IC, HOPEI i let n susceptibles and a single 
infective be there in a population of size n+1. let p(r,t) 
denote the probability of there being r susceptibles in the 
system at any time t. Assuming i infections in (t,t + At) 

(i s= 1,2,... ,m) with p 15 pg,..., P m the infection rates, let 
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p^r(n+l~r) At, -f3gr (n+l-r) At, . .., p m r(n+l-r) At be the 
probabilities of one infection, two infections, three infections, 
m infections when the number of susceptibles is r at any 
time. When we write down the relation between probabilities of 
different states of the system, we are unable to get the 
differential- difference equations system as compact as in 
the case of simple epidemic model with single infection. Here 
the system is very much dependent on m and n. Of course we 
assume that m < n with m very small in comparison to n. 


jO *j 

lor the case m <_ n/2 when n is even or m <_ when 

n is odd, we get a differential-difference equation system 

in which the first m equations (r = 0, 1, 2, . . . ,m-l) and the 

last m equations (r = n-m+l,...,n) are written separately with 

great care on the relations between probabilities of different 

states. But the middle n+1 - 2m equations can be represented 

by a single equation with r varying from m to n-m. Changing 

Pi . 

the time-scale to x = t(3^ and writing G ^ = p— (1= 2,3, 
we get the following system of differential-difference equation 


SLOjJI — np(l, t) + Gg 2(n-l)p (2, f)+. • »+m(n— hh-1)g^ p(m, t ) 

MAjJEI _ 2(n-l)p (2, r)+3(n-2)Gg p(3, t)+. ..+(m+l)(n-m)G m p(m+l, t) 

- n p(l, t). 


= m(n+l-m)p (m, T)+(m+l) (n-m)Gg p(m+l, T ) + ... 

+ ( 2m- 1 ) (n+2-2m )G m p ( 2m- 1, t ) - ( 1+G 2+S3+ * • * +e ni-l ^ 

(m-l)(n-m+2) p(m-l,x) 



= ( r +l)(n-r)p(r+l, T)+(r+2)(n-r_l)e 2 p(r+2, t) 

+ (r+3)(n-r_2)e 3 p (r+3, t)+. . .+(r+m)(n~r-m+l)e m p(r+m, t) 
_(l+e 2 +e 3 +. . .+e m ) r(n+l-r) p(r,r) 

r = . . . ,n-m 

= (n-m+2 ) (m-1 ) p (n_m+2 , r)+e 2 (n-m+3 ) (m-2 )p (n-m+3 , t) 

+...+G _ 1 np(n,T)-(l+e 2 +8 3 +...+e m )(n-m+l)mp(ia«m+l,T) 


itoiSrtiijJ. _ np (n, t) - ( f+Gp+Gg-t-. • *+€ m ) 2(n-l) p(n-l,r) 

^piPxlL = _(i+e 2 + e 3 + ... + e m ) np(n,T) 

with initial conditions 

p(r,0) = 1 when r = n 


( 8 ) 


(9) 


n-1 


p(r,0) = t) when r ^ n. 

Por the case m > n/2 when n is even or m > - g— when 
n is odd, we have to write down the n+1 equations separately. 
In this case no two equations can be combined into single 
equation. Following the previous procedure, we get the 

differentiated if ference equation system as follows. 


_ np(l, r)+2(n-l)e 2 p(2, t) +...+ m(n-m+l) e m P( m >0 


SPlP-jLil* 2(n-l)p(2, x)+3(n-2)e 2 p(3, t) + . . . + 6 m (m+l)(n_m)p (b*1, t) 

d*t v 

- np(l.r) 



112 


dp (riul, t) / , w 6g P(m+ 1 5 T ) +••• 

" dT ^ = m (n+l-m) p(m, t)+C®+ a oW , . 

+ e n-m + l n P^ T) " C 2 

AnZ-rn ^ s < rtu.2 ) ^ 6 2 P(®+2,t)+... 

dV^ = ( m+ l)( n ~ m ) p(bi+ 1 j t ^ + ^ 

„ .p + ...+e m ) m ( n -m+l) p(m,T) 

+ e n-m n P( n 5' r )"( 1+e 2 +e 3 

win N ^/ ra+ 3)(n-o-2)e ? p(m+3, t)+... 

_ ( m+ £ ) (n-iTU l)p (ff -+ 2 5 

. v ri+ e 9 +e.+ --' +G m )(m+1)(n “ m)p(mfljT) 
+e n-^l n P (n ’ T)_( ^ 2 5 


(9) 


&igri-ll = np(n, t ) - (l+e g +e 3 + ' 


,p ) 2(n-l) p(n-l,T) 
• ♦ ♦ +b -m/ 


= _ ( i + 6p + e2 + ... + e m 


ej »p( n ’ T ) 


with initial conditions p( r J G -^ 


1 when r = n otherwise o 


(10) 


5.2.3 !EWO SPECIAL CASES * 


(i) n = 3, m = 2 i* e * m 


s n-1 

< 11 ^ ^ 


MPjjl 

dx 

dr 


3p(l, x) + 4G g pt 2 *"^ 

4p(2,x) + 30 p p( 3 > t) " ^ 


( 11 ) 


. d P.(.2jjl _ 3p (3, t) - 4(1+Gg) P( 2 » T ' ) 

— —3 ( 1+6 y ) p ( 3 j x ) 

p( 3 ,o) = 1 , P(r.°) = ° 1 * 


(12) 
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(ii) n = 6, m = 2 i.e. m < n/2 
_ 6p(l,x) + loe g p (2, t) 

= lOp(2,x) + 12Sg p(3,x) - 6p(l,x) 

M|xJ 1 = I2p(3,x) + 12Gg p(4,t) - I0(l+e 2 ) p(2,x) 

= 12p(4,x) + loe 2 p(5,x) .. 12(l+6g) p(3,x) 

_ I0p(5, x) + 6Gg p(6, t) - 12(l+Gg) p(4,x) 

MJkll = 6p(6, x) - 10 ( 1+G g ) p ( 5, t) 
dp>X|-sJL}. _ „6(1+Gg) p(6,x) 

p(6,0) = 1 otherwise 0. (14) 

for any value of n and m given, v/e a get a system of 
differential-difference equations of first order. Inis system, 
with initial conditions can be solved numerically by "ill Method. 
Applying Gill-Method, we solve the two special cases numerically. 

5.2.4 GRAPHS 

for n = 6, m = 2 

(i) In figures 5.3 - 5.9, we have plotted p(0, t), p(l,x), 
p(2,x), p (3 , x) , p (4, x) , p(5,x) and p (6, t) respectively 
against time for Gg = 0, 0*2 and 0.5. 

(ii) In Pig. 5.10, we have plotted all probabilities against 

time for Gg = 0. 
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(iii) In Fig. 5.11, we have plotted all probabilities 
against time for e 2 = 0.2. 

(iv) In Fig. 5.12, we have plotted all probabilities 
against time for e 2 = 0.5. 

(v) In Fig. 5.13, we have plotted mean against time for 
S 2 — 0, 0.2 an u 0.0# 

(vi) In Fig. 5.14, we have plotted variance against time 
for 6g = 0, 0.2 and 0.5. 

For n = 3, m = 2 

(i) In Figures 5.15 - 5.18, we have plotted p (0, r), p(l,r), 
p(2,r) and p(3,t) respectively against time for e 2 = 0, 
0.2 and 0.5. 

(ii) In Fig. 5.19, we nave plotted all probabilities 
against time for e 2 = 0. 

(iii) In Fig. 5 .20, we have plotted all probabilities against 
time for 0 2 = 0.2. 

( iv) In Fig. 5.21, we have plotted all probabilities against 
time for e 2 = 0.5. 

(v) In Fig. 5.22, we have plotted mean against time for 
Og = 0, 0.2 and 0.5. 

(vi) In Fig. 5.23, we have plotted variance against time 

$ 

for e Q = 0, 0.2 and 0.5. 

We observe, as expected that multiple infection causes 
the process of infection to be faster and it takes little time 
for all the susceptibles to be infected. The probabilities for 
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which r is near to n are being zero faster and for which r is 
close to zero are reaching 1 fastly. 


5.3 G1EERA1 EPI1ELIIC H0DE1 


5.3.1 ST OGiiA S TIG MODEL 

In the population size of n+a where n and a denote the 
initial number of susceptibles and infectives respectively, we 
assume i infections (i = l,2,...,m) and single removal in the 
time interval (t, t + At) with p^Pg,...^ as the infection 
rates and Y the removal rate, let p (t) denote the probability 

X! 3 

of there being r susceptibles and s infectives in the system 

at any time t. let p^_ rs &t (i = l,2,...,m) denote the 

probability of i infections and ys At, the probability of 

single removal in time interval (t, t + At). Changing the time— 

Y P -v 

scale to t = p-jt and writing = P and Gj_ = pp (i=2 ,3, , -m V 
we get the following system of differential-differential 




equations 

d P rs ( T ) 
“*'dx ~ 


(r+l)(s-l) P r+1;S _i( T ) + e s (r+Z)(s-Z) P r+ 2 t s_gC ) 


+ ...-t€ m (r+m)(s-m)p r+m) s_ In (T) + P r)SH .i<h “(s+l) 

- {(i+e 9 +e*+...+ e )rs + ps> p rq ( T ) 

25 min (r,m) rs 

r = 0,1,2, , ..,n-m 

s — m, m+ 1 , . . * , m+ ( a— 1 ) 


wi 


ith initial conditions 


n (o) = 1, otherwise 0 
•^na ' 1 


( 2 ) 
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and 0 < r < n, 0 < s < n+a, 0 < r + s < n+a (3) 

With the above restrictions (3) on r, s and r+s, we 

get (n+l)(n+a+l) - number of ncn-zero probabilities. 

Of y/hich a(n~m+l) number of different ial- difference equations 
can be obtained directly from the above equation (1), 
substituting restricted values of r and s. 

Even all other probabilities can be obtained from (1) 
if we make the following convention. 

(i) Probability with suffix negative will be taken as zero. 

(ii) Condition (3) will be taken care of. 

5.3.2 SPECIA L CA SE : 

We present here as an example a complete system of 
differential- difference equations for n = 3, a = 2 and m = 2. 

In this case we get a system of 18 differential-difference 
equations as follows. 


fjoo 

~dV 


PP 


ol 


?fol 
~*dT 


= -P P 0 i + 2 P P 


o2 


ffog. 

“Sir*"' 


“ 2p -P 0 2 + P 11 + 3P P o3 


= “ 3p p o3 + 2Pp l2 + 2e 2 P 21 + P o4 
-T?’ = “ 4p p o4 + 3p 13 + 4 ®2 p 22 + 5C P o5 
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4p n 3 

*Tx~“ = “ 5p p o5 + 4p l4 + 6G 2 p 23 

p Pll 

-(1+ p ) Pn + P ±2 
-(2+2P) p^2 + G P2l + P 13 

■4^* = -( 3 + 3p ) P 13 + 4p 22 + 3e 2 P 31 + 4p P 14 
dp 1/L 

~ctr = " (4 + 4p) p 14 + 6p 23 + 6G 2 p 32 
pp 2 1 

-(2 + 2Gg + p) P24 + G p P22 
-(4 + 4e 2 + 2 P) p 2 g + 3 P3p + 3p P 23 

dp o rz v r 

■'“d'r" = ~(6 + 66 2 + 3p) P23 + 6p 32 

dl?30 n n 

“Tt“ = p P31 

r\v’\ 

=~( 3 + 3 e + p) p 3 i + 2p P32 

(It) r - Q 

= -(6 + 6 e + 2p ) P 32 

We solve this system of first order differential- difference 
equations nuner ically by SiU-Hethod for different values of 

and e 2 * 


5f20 = 

"dr 
dp 2 i _ 

~'cTt " “ 

dp 22 
"dV ~ 


4p -j_Q 

'"*d’V" = 

11 

"d'x 

dp l2 

~aV’ ~ 
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5 .3 ,3 G-lLAPi-IS 

(i) In Pig. 5.24, we have plotted P 00 ( t), the probability 
that all suscoyt ibles get infected and all infect ives 
are removed. 

(ii) In Pig. 5.25 - 5.29, we have plotted all probabilities 
of extinction for suscep titles, 

(iii) In Pig. 5.30, we have plotted p 32 (t), ’ fc5ie P^^ability 
that no susceptible gets infected and no infective 

is removed. 
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Pq2 ^ Sjj = 0*5 p = 0.5 

2 e 2 = 0.5 p = 0.1 

3 e 2 = 0,1 P = 0.5 

4 0.5 P = o.l 


0.5 - 
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